|
12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307130813091310131113121313131413151316131713181319132013211322132313241325132613271328132913301331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400140114021403140414051406140714081409141014111412141314141415141614171418141914201421142214231424142514261427142814291430143114321433143414351436143714381439144014411442144314441445144614471448144914501451145214531454145514561457145814591460146114621463146414651466146714681469147014711472147314741475147614771478147914801481148214831484148514861487148814891490149114921493149414951496149714981499150015011502150315041505150615071508150915101511151215131514151515161517151815191520152115221523152415251526152715281529153015311532153315341535153615371538153915401541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610161116121613161416151616161716181619162016211622162316241625162616271628162916301631163216331634163516361637163816391640164116421643164416451646164716481649165016511652165316541655165616571658165916601661166216631664166516661667166816691670167116721673167416751676167716781679168016811682168316841685168616871688168916901691169216931694169516961697169816991700170117021703170417051706170717081709171017111712171317141715171617171718171917201721172217231724172517261727172817291730173117321733173417351736173717381739174017411742174317441745174617471748174917501751175217531754175517561757175817591760176117621763176417651766176717681769177017711772177317741775177617771778177917801781178217831784178517861787178817891790179117921793179417951796179717981799180018011802180318041805180618071808180918101811181218131814181518161817181818191820182118221823182418251826182718281829183018311832183318341835183618371838183918401841184218431844184518461847184818491850185118521853185418551856185718581859186018611862186318641865186618671868186918701871187218731874187518761877187818791880188118821883188418851886188718881889189018911892189318941895189618971898189919001901190219031904190519061907190819091910191119121913191419151916191719181919192019211922192319241925192619271928192919301931193219331934193519361937193819391940194119421943194419451946194719481949195019511952195319541955195619571958195919601961196219631964196519661967196819691970197119721973197419751976197719781979198019811982198319841985198619871988198919901991199219931994199519961997199819992000200120022003200420052006200720082009201020112012201320142015201620172018201920202021202220232024202520262027202820292030203120322033203420352036203720382039204020412042204320442045204620472048204920502051205220532054205520562057205820592060206120622063206420652066206720682069207020712072207320742075207620772078207920802081208220832084208520862087208820892090209120922093209420952096209720982099210021012102210321042105210621072108210921102111211221132114211521162117211821192120212121222123212421252126212721282129213021312132213321342135213621372138213921402141214221432144214521462147214821492150215121522153215421552156215721582159216021612162216321642165216621672168216921702171217221732174217521762177217821792180218121822183218421852186218721882189219021912192219321942195219621972198219922002201220222032204220522062207220822092210221122122213221422152216221722182219222022212222222322242225222622272228222922302231223222332234223522362237223822392240224122422243224422452246224722482249225022512252225322542255225622572258225922602261226222632264226522662267226822692270227122722273227422752276227722782279228022812282228322842285228622872288228922902291229222932294229522962297229822992300230123022303230423052306230723082309231023112312231323142315231623172318231923202321232223232324232523262327232823292330233123322333233423352336233723382339234023412342234323442345234623472348234923502351235223532354235523562357235823592360236123622363236423652366236723682369237023712372237323742375237623772378237923802381238223832384238523862387238823892390239123922393239423952396239723982399240024012402240324042405240624072408240924102411241224132414241524162417241824192420242124222423242424252426242724282429243024312432243324342435243624372438243924402441244224432444244524462447244824492450245124522453245424552456245724582459246024612462246324642465246624672468246924702471247224732474247524762477247824792480248124822483248424852486248724882489249024912492249324942495249624972498249925002501250225032504250525062507250825092510251125122513251425152516251725182519252025212522252325242525252625272528252925302531253225332534253525362537253825392540254125422543254425452546254725482549255025512552255325542555255625572558255925602561256225632564256525662567256825692570257125722573257425752576257725782579258025812582258325842585258625872588258925902591259225932594259525962597259825992600260126022603260426052606260726082609261026112612261326142615261626172618261926202621262226232624262526262627262826292630263126322633263426352636263726382639264026412642264326442645264626472648264926502651265226532654265526562657265826592660266126622663266426652666266726682669267026712672267326742675267626772678267926802681268226832684268526862687268826892690269126922693269426952696269726982699270027012702270327042705270627072708270927102711271227132714271527162717271827192720272127222723272427252726272727282729273027312732273327342735273627372738273927402741274227432744274527462747274827492750275127522753275427552756275727582759276027612762276327642765276627672768276927702771277227732774277527762777277827792780278127822783278427852786278727882789279027912792279327942795279627972798279928002801280228032804280528062807280828092810281128122813281428152816281728182819282028212822282328242825282628272828282928302831283228332834283528362837283828392840284128422843284428452846284728482849285028512852285328542855285628572858285928602861286228632864286528662867286828692870287128722873287428752876287728782879288028812882288328842885288628872888288928902891289228932894289528962897289828992900290129022903290429052906290729082909291029112912291329142915291629172918291929202921292229232924292529262927292829292930293129322933293429352936293729382939294029412942294329442945294629472948294929502951295229532954295529562957295829592960296129622963296429652966296729682969297029712972297329742975297629772978297929802981298229832984298529862987298829892990299129922993299429952996299729982999300030013002300330043005300630073008300930103011301230133014301530163017301830193020302130223023302430253026302730283029303030313032303330343035303630373038303930403041304230433044304530463047304830493050305130523053305430553056305730583059306030613062306330643065306630673068306930703071307230733074307530763077307830793080308130823083308430853086308730883089309030913092309330943095309630973098309931003101310231033104310531063107310831093110311131123113311431153116311731183119312031213122312331243125312631273128312931303131313231333134313531363137313831393140314131423143314431453146314731483149315031513152315331543155315631573158315931603161316231633164316531663167316831693170317131723173317431753176317731783179318031813182318331843185318631873188318931903191319231933194319531963197319831993200320132023203320432053206320732083209321032113212321332143215321632173218321932203221322232233224322532263227322832293230323132323233323432353236323732383239324032413242324332443245324632473248324932503251325232533254325532563257325832593260326132623263326432653266326732683269327032713272327332743275327632773278327932803281328232833284328532863287328832893290329132923293329432953296329732983299330033013302330333043305330633073308330933103311331233133314331533163317331833193320332133223323332433253326332733283329333033313332333333343335333633373338333933403341334233433344334533463347334833493350335133523353335433553356335733583359336033613362336333643365336633673368336933703371337233733374337533763377337833793380338133823383338433853386338733883389339033913392339333943395339633973398339934003401340234033404340534063407340834093410341134123413341434153416341734183419342034213422342334243425342634273428342934303431343234333434343534363437343834393440344134423443344434453446344734483449345034513452345334543455345634573458345934603461346234633464346534663467346834693470347134723473347434753476347734783479348034813482348334843485348634873488348934903491349234933494349534963497349834993500350135023503350435053506350735083509351035113512351335143515351635173518351935203521352235233524352535263527352835293530353135323533353435353536353735383539354035413542354335443545354635473548354935503551355235533554355535563557355835593560356135623563356435653566356735683569357035713572357335743575357635773578357935803581358235833584358535863587358835893590359135923593359435953596359735983599360036013602360336043605360636073608360936103611361236133614361536163617361836193620362136223623362436253626362736283629363036313632363336343635363636373638363936403641364236433644364536463647364836493650365136523653365436553656 |
- #include <math.h>
- #include <stdlib.h>
- #include <string.h>
- #include <stdio.h>
- #include <complex.h>
- #ifdef complex
- #undef complex
- #endif
- #ifdef I
- #undef I
- #endif
-
- #include "common.h"
-
- typedef blasint integer;
-
- typedef unsigned int uinteger;
- typedef char *address;
- typedef short int shortint;
- typedef float real;
- typedef double doublereal;
- typedef struct { real r, i; } complex;
- typedef struct { doublereal r, i; } doublecomplex;
- #ifdef _MSC_VER
- static inline _Fcomplex Cf(complex *z) {_Fcomplex zz={z->r , z->i}; return zz;}
- static inline _Dcomplex Cd(doublecomplex *z) {_Dcomplex zz={z->r , z->i};return zz;}
- static inline _Fcomplex * _pCf(complex *z) {return (_Fcomplex*)z;}
- static inline _Dcomplex * _pCd(doublecomplex *z) {return (_Dcomplex*)z;}
- #else
- static inline _Complex float Cf(complex *z) {return z->r + z->i*_Complex_I;}
- static inline _Complex double Cd(doublecomplex *z) {return z->r + z->i*_Complex_I;}
- static inline _Complex float * _pCf(complex *z) {return (_Complex float*)z;}
- static inline _Complex double * _pCd(doublecomplex *z) {return (_Complex double*)z;}
- #endif
- #define pCf(z) (*_pCf(z))
- #define pCd(z) (*_pCd(z))
- typedef int logical;
- typedef short int shortlogical;
- typedef char logical1;
- typedef char integer1;
-
- #define TRUE_ (1)
- #define FALSE_ (0)
-
- /* Extern is for use with -E */
- #ifndef Extern
- #define Extern extern
- #endif
-
- /* I/O stuff */
-
- typedef int flag;
- typedef int ftnlen;
- typedef int ftnint;
-
- /*external read, write*/
- typedef struct
- { flag cierr;
- ftnint ciunit;
- flag ciend;
- char *cifmt;
- ftnint cirec;
- } cilist;
-
- /*internal read, write*/
- typedef struct
- { flag icierr;
- char *iciunit;
- flag iciend;
- char *icifmt;
- ftnint icirlen;
- ftnint icirnum;
- } icilist;
-
- /*open*/
- typedef struct
- { flag oerr;
- ftnint ounit;
- char *ofnm;
- ftnlen ofnmlen;
- char *osta;
- char *oacc;
- char *ofm;
- ftnint orl;
- char *oblnk;
- } olist;
-
- /*close*/
- typedef struct
- { flag cerr;
- ftnint cunit;
- char *csta;
- } cllist;
-
- /*rewind, backspace, endfile*/
- typedef struct
- { flag aerr;
- ftnint aunit;
- } alist;
-
- /* inquire */
- typedef struct
- { flag inerr;
- ftnint inunit;
- char *infile;
- ftnlen infilen;
- ftnint *inex; /*parameters in standard's order*/
- ftnint *inopen;
- ftnint *innum;
- ftnint *innamed;
- char *inname;
- ftnlen innamlen;
- char *inacc;
- ftnlen inacclen;
- char *inseq;
- ftnlen inseqlen;
- char *indir;
- ftnlen indirlen;
- char *infmt;
- ftnlen infmtlen;
- char *inform;
- ftnint informlen;
- char *inunf;
- ftnlen inunflen;
- ftnint *inrecl;
- ftnint *innrec;
- char *inblank;
- ftnlen inblanklen;
- } inlist;
-
- #define VOID void
-
- union Multitype { /* for multiple entry points */
- integer1 g;
- shortint h;
- integer i;
- /* longint j; */
- real r;
- doublereal d;
- complex c;
- doublecomplex z;
- };
-
- typedef union Multitype Multitype;
-
- struct Vardesc { /* for Namelist */
- char *name;
- char *addr;
- ftnlen *dims;
- int type;
- };
- typedef struct Vardesc Vardesc;
-
- struct Namelist {
- char *name;
- Vardesc **vars;
- int nvars;
- };
- typedef struct Namelist Namelist;
-
- #define abs(x) ((x) >= 0 ? (x) : -(x))
- #define dabs(x) (fabs(x))
- #define f2cmin(a,b) ((a) <= (b) ? (a) : (b))
- #define f2cmax(a,b) ((a) >= (b) ? (a) : (b))
- #define dmin(a,b) (f2cmin(a,b))
- #define dmax(a,b) (f2cmax(a,b))
- #define bit_test(a,b) ((a) >> (b) & 1)
- #define bit_clear(a,b) ((a) & ~((uinteger)1 << (b)))
- #define bit_set(a,b) ((a) | ((uinteger)1 << (b)))
-
- #define abort_() { sig_die("Fortran abort routine called", 1); }
- #define c_abs(z) (cabsf(Cf(z)))
- #define c_cos(R,Z) { pCf(R)=ccos(Cf(Z)); }
- #ifdef _MSC_VER
- #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]);}
- #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]);}
- #else
- #define c_div(c, a, b) {pCf(c) = Cf(a)/Cf(b);}
- #define z_div(c, a, b) {pCd(c) = Cd(a)/Cd(b);}
- #endif
- #define c_exp(R, Z) {pCf(R) = cexpf(Cf(Z));}
- #define c_log(R, Z) {pCf(R) = clogf(Cf(Z));}
- #define c_sin(R, Z) {pCf(R) = csinf(Cf(Z));}
- //#define c_sqrt(R, Z) {*(R) = csqrtf(Cf(Z));}
- #define c_sqrt(R, Z) {pCf(R) = csqrtf(Cf(Z));}
- #define d_abs(x) (fabs(*(x)))
- #define d_acos(x) (acos(*(x)))
- #define d_asin(x) (asin(*(x)))
- #define d_atan(x) (atan(*(x)))
- #define d_atn2(x, y) (atan2(*(x),*(y)))
- #define d_cnjg(R, Z) { pCd(R) = conj(Cd(Z)); }
- #define r_cnjg(R, Z) { pCf(R) = conjf(Cf(Z)); }
- #define d_cos(x) (cos(*(x)))
- #define d_cosh(x) (cosh(*(x)))
- #define d_dim(__a, __b) ( *(__a) > *(__b) ? *(__a) - *(__b) : 0.0 )
- #define d_exp(x) (exp(*(x)))
- #define d_imag(z) (cimag(Cd(z)))
- #define r_imag(z) (cimagf(Cf(z)))
- #define d_int(__x) (*(__x)>0 ? floor(*(__x)) : -floor(- *(__x)))
- #define r_int(__x) (*(__x)>0 ? floor(*(__x)) : -floor(- *(__x)))
- #define d_lg10(x) ( 0.43429448190325182765 * log(*(x)) )
- #define r_lg10(x) ( 0.43429448190325182765 * log(*(x)) )
- #define d_log(x) (log(*(x)))
- #define d_mod(x, y) (fmod(*(x), *(y)))
- #define u_nint(__x) ((__x)>=0 ? floor((__x) + .5) : -floor(.5 - (__x)))
- #define d_nint(x) u_nint(*(x))
- #define u_sign(__a,__b) ((__b) >= 0 ? ((__a) >= 0 ? (__a) : -(__a)) : -((__a) >= 0 ? (__a) : -(__a)))
- #define d_sign(a,b) u_sign(*(a),*(b))
- #define r_sign(a,b) u_sign(*(a),*(b))
- #define d_sin(x) (sin(*(x)))
- #define d_sinh(x) (sinh(*(x)))
- #define d_sqrt(x) (sqrt(*(x)))
- #define d_tan(x) (tan(*(x)))
- #define d_tanh(x) (tanh(*(x)))
- #define i_abs(x) abs(*(x))
- #define i_dnnt(x) ((integer)u_nint(*(x)))
- #define i_len(s, n) (n)
- #define i_nint(x) ((integer)u_nint(*(x)))
- #define i_sign(a,b) ((integer)u_sign((integer)*(a),(integer)*(b)))
- #define pow_dd(ap, bp) ( pow(*(ap), *(bp)))
- #define pow_si(B,E) spow_ui(*(B),*(E))
- #define pow_ri(B,E) spow_ui(*(B),*(E))
- #define pow_di(B,E) dpow_ui(*(B),*(E))
- #define pow_zi(p, a, b) {pCd(p) = zpow_ui(Cd(a), *(b));}
- #define pow_ci(p, a, b) {pCf(p) = cpow_ui(Cf(a), *(b));}
- #define pow_zz(R,A,B) {pCd(R) = cpow(Cd(A),*(B));}
- #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++ = ' '; }
- #define s_cmp(a,b,c,d) ((integer)strncmp((a),(b),f2cmin((c),(d))))
- #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]; }
- #define sig_die(s, kill) { exit(1); }
- #define s_stop(s, n) {exit(0);}
- #define z_abs(z) (cabs(Cd(z)))
- #define z_exp(R, Z) {pCd(R) = cexp(Cd(Z));}
- #define z_sqrt(R, Z) {pCd(R) = csqrt(Cd(Z));}
- #define myexit_() break;
- #define mycycle_() continue;
- #define myceiling_(w) {ceil(w)}
- #define myhuge_(w) {HUGE_VAL}
- //#define mymaxloc_(w,s,e,n) {if (sizeof(*(w)) == sizeof(double)) dmaxloc_((w),*(s),*(e),n); else dmaxloc_((w),*(s),*(e),n);}
- #define mymaxloc_(w,s,e,n) dmaxloc_(w,*(s),*(e),n)
-
- /* procedure parameter types for -A and -C++ */
-
- #define F2C_proc_par_types 1
- #ifdef __cplusplus
- typedef logical (*L_fp)(...);
- #else
- typedef logical (*L_fp)();
- #endif
- #if 0
- static float spow_ui(float x, integer n) {
- float pow=1.0; unsigned long int u;
- if(n != 0) {
- if(n < 0) n = -n, x = 1/x;
- for(u = n; ; ) {
- if(u & 01) pow *= x;
- if(u >>= 1) x *= x;
- else break;
- }
- }
- return pow;
- }
- static double dpow_ui(double x, integer n) {
- double pow=1.0; unsigned long int u;
- if(n != 0) {
- if(n < 0) n = -n, x = 1/x;
- for(u = n; ; ) {
- if(u & 01) pow *= x;
- if(u >>= 1) x *= x;
- else break;
- }
- }
- return pow;
- }
- #ifdef _MSC_VER
- static _Fcomplex cpow_ui(complex x, integer n) {
- complex pow={1.0,0.0}; unsigned long int u;
- if(n != 0) {
- if(n < 0) n = -n, x.r = 1/x.r, x.i=1/x.i;
- for(u = n; ; ) {
- if(u & 01) pow.r *= x.r, pow.i *= x.i;
- if(u >>= 1) x.r *= x.r, x.i *= x.i;
- else break;
- }
- }
- _Fcomplex p={pow.r, pow.i};
- return p;
- }
- #else
- static _Complex float cpow_ui(_Complex float x, integer n) {
- _Complex float pow=1.0; unsigned long int u;
- if(n != 0) {
- if(n < 0) n = -n, x = 1/x;
- for(u = n; ; ) {
- if(u & 01) pow *= x;
- if(u >>= 1) x *= x;
- else break;
- }
- }
- return pow;
- }
- #endif
- #ifdef _MSC_VER
- static _Dcomplex zpow_ui(_Dcomplex x, integer n) {
- _Dcomplex pow={1.0,0.0}; unsigned long int u;
- if(n != 0) {
- if(n < 0) n = -n, x._Val[0] = 1/x._Val[0], x._Val[1] =1/x._Val[1];
- for(u = n; ; ) {
- if(u & 01) pow._Val[0] *= x._Val[0], pow._Val[1] *= x._Val[1];
- if(u >>= 1) x._Val[0] *= x._Val[0], x._Val[1] *= x._Val[1];
- else break;
- }
- }
- _Dcomplex p = {pow._Val[0], pow._Val[1]};
- return p;
- }
- #else
- static _Complex double zpow_ui(_Complex double x, integer n) {
- _Complex double pow=1.0; unsigned long int u;
- if(n != 0) {
- if(n < 0) n = -n, x = 1/x;
- for(u = n; ; ) {
- if(u & 01) pow *= x;
- if(u >>= 1) x *= x;
- else break;
- }
- }
- return pow;
- }
- #endif
- static integer pow_ii(integer x, integer n) {
- integer pow; unsigned long int u;
- if (n <= 0) {
- if (n == 0 || x == 1) pow = 1;
- else if (x != -1) pow = x == 0 ? 1/x : 0;
- else n = -n;
- }
- if ((n > 0) || !(n == 0 || x == 1 || x != -1)) {
- u = n;
- for(pow = 1; ; ) {
- if(u & 01) pow *= x;
- if(u >>= 1) x *= x;
- else break;
- }
- }
- return pow;
- }
- static integer dmaxloc_(double *w, integer s, integer e, integer *n)
- {
- double m; integer i, mi;
- for(m=w[s-1], mi=s, i=s+1; i<=e; i++)
- if (w[i-1]>m) mi=i ,m=w[i-1];
- return mi-s+1;
- }
- static integer smaxloc_(float *w, integer s, integer e, integer *n)
- {
- float m; integer i, mi;
- for(m=w[s-1], mi=s, i=s+1; i<=e; i++)
- if (w[i-1]>m) mi=i ,m=w[i-1];
- return mi-s+1;
- }
- #endif
- /* -- translated by f2c (version 20000121).
- You must link the resulting object file with the libraries:
- -lf2c -lm (in that order)
- */
-
-
-
- /* Common Block Declarations */
-
- struct {
- integer infot, noutc;
- logical ok;
- } infoc_;
-
- #define infoc_1 infoc_
-
- struct {
- char srnamt[12];
- } srnamc_;
-
- #define srnamc_1 srnamc_
-
- /* Table of constant values */
-
- static integer c__1 = 1;
- static integer c__65 = 65;
- static doublereal c_b90 = 1.;
- static doublereal c_b104 = 0.;
- static integer c__6 = 6;
- static logical c_true = TRUE_;
- static integer c__0 = 0;
- static logical c_false = FALSE_;
-
- /* Main program MAIN__() */ int main()
- {
- /* Initialized data */
-
- static char snames[6][13] = {"cblas_dgemm ", "cblas_dsymm ", "cblas_dtrmm ", "cblas_dtrsm ", "cblas_dsyrk ", "cblas_dsyr2k"};
-
- /* System generated locals */
- integer i__1, i__2, i__3;
- doublereal d__1;
-
- /* Builtin functions */
- integer s_rsle(), do_lio(), e_rsle(), f_open(), s_wsfe(), do_fio(),
- e_wsfe(), s_wsle(), e_wsle(), s_rsfe(), e_rsfe();
- integer f_clos();
-
- /* Local variables */
- static integer nalf, idim[9];
- static logical same;
- static integer nbet, ntra;
- static logical rewi;
- extern /* Subroutine */ int dchk1_(), dchk2_(), dchk3_(), dchk4_(),
- dchk5_();
- static doublereal c__[4225] /* was [65][65] */, g[65];
- static integer i__, j;
- extern doublereal ddiff_();
- static integer n;
- static logical fatal;
- static doublereal w[130];
- extern /* Subroutine */ int dmmch_();
- static logical trace;
- static integer nidim;
- static char snaps[32];
- static integer isnum;
- static logical ltest[6];
- static doublereal aa[4225], ab[8450] /* was [65][130] */, bb[4225],
- cc[4225], as[4225], bs[4225], cs[4225], ct[65];
- static logical sfatal, corder;
- static char snamet[12], transa[1], transb[1];
- static doublereal thresh;
- static logical rorder;
- extern /* Subroutine */ int cd3chke_();
- static integer layout;
- static logical ltestt, tsterr;
- static doublereal alf[7];
- extern logical lde_();
- static doublereal bet[7], eps, err;
- char tmpchar;
-
- /* Test program for the DOUBLE PRECISION Level 3 Blas. */
-
- /* The program must be driven by a short data file. The first 13 records */
- /* of the file are read using list-directed input, the last 6 records */
- /* are read using the format ( A12, L2 ). An annotated example of a data */
- /* file can be obtained by deleting the first 3 characters from the */
- /* following 19 lines: */
- /* 'DBLAT3.SNAP' NAME OF SNAPSHOT OUTPUT FILE */
- /* -1 UNIT NUMBER OF SNAPSHOT FILE (NOT USED IF .LT. 0) */
- /* F LOGICAL FLAG, T TO REWIND SNAPSHOT FILE AFTER EACH RECORD. */
- /* F LOGICAL FLAG, T TO STOP ON FAILURES. */
- /* T LOGICAL FLAG, T TO TEST ERROR EXITS. */
- /* 2 0 TO TEST COLUMN-MAJOR, 1 TO TEST ROW-MAJOR, 2 TO TEST BOTH */
- /* 16.0 THRESHOLD VALUE OF TEST RATIO */
- /* 6 NUMBER OF VALUES OF N */
- /* 0 1 2 3 5 9 VALUES OF N */
- /* 3 NUMBER OF VALUES OF ALPHA */
- /* 0.0 1.0 0.7 VALUES OF ALPHA */
- /* 3 NUMBER OF VALUES OF BETA */
- /* 0.0 1.0 1.3 VALUES OF BETA */
- /* cblas_dgemm T PUT F FOR NO TEST. SAME COLUMNS. */
- /* cblas_dsymm T PUT F FOR NO TEST. SAME COLUMNS. */
- /* cblas_dtrmm T PUT F FOR NO TEST. SAME COLUMNS. */
- /* cblas_dtrsm T PUT F FOR NO TEST. SAME COLUMNS. */
- /* cblas_dsyrk T PUT F FOR NO TEST. SAME COLUMNS. */
- /* cblas_dsyr2k T PUT F FOR NO TEST. SAME COLUMNS. */
-
- /* See: */
-
- /* Dongarra J. J., Du Croz J. J., Duff I. S. and Hammarling S. */
- /* A Set of Level 3 Basic Linear Algebra Subprograms. */
-
- /* Technical Memorandum No.88 (Revision 1), Mathematics and */
- /* Computer Science Division, Argonne National Laboratory, 9700 */
- /* South Cass Avenue, Argonne, Illinois 60439, US. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Local Scalars .. */
- /* .. Local Arrays .. */
- /* .. External Functions .. */
- /* .. External Subroutines .. */
- /* .. Intrinsic Functions .. */
- /* .. Scalars in Common .. */
- /* .. Common blocks .. */
- /* .. Data statements .. */
- /* .. Executable Statements .. */
-
- /* Read name and unit number for summary output file and open file. */
-
- infoc_1.noutc = 6;
- /* Read name and unit number for snapshot output file and open file. */
-
- char line[80];
-
- fgets(line,80,stdin);
- sscanf(line,"'%s'",snaps);
- fgets(line,80,stdin);
- #ifdef USE64BITINT
- sscanf(line,"%ld",&ntra);
- #else
- sscanf(line,"%d",&ntra);
- #endif
- trace = ntra >= 0;
- if (trace) {
- /* o__1.oerr = 0;
- o__1.ounit = ntra;
- o__1.ofnmlen = 32;
- o__1.ofnm = snaps;
- o__1.orl = 0;
- o__1.osta = "NEW";
- o__1.oacc = 0;
- o__1.ofm = 0;
- o__1.oblnk = 0;
- f_open(&o__1);*/
- }
- /* Read the flag that directs rewinding of the snapshot file. */
- fgets(line,80,stdin);
- sscanf(line,"%d",&rewi);
- rewi = rewi && trace;
- /* Read the flag that directs stopping on any failure. */
- fgets(line,80,stdin);
- sscanf(line,"%c",&tmpchar);
- /* Read the flag that indicates whether error exits are to be tested. */
- sfatal=FALSE_;
- if (tmpchar=='T')sfatal=TRUE_;
- fgets(line,80,stdin);
- sscanf(line,"%c",&tmpchar);
- /* Read the flag that indicates whether error exits are to be tested. */
- tsterr=FALSE_;
- if (tmpchar=='T')tsterr=TRUE_;
- /* Read the flag that indicates whether row-major data layout to be tested. */
- fgets(line,80,stdin);
- sscanf(line,"%d",&layout);
- /* Read the threshold value of the test ratio */
- fgets(line,80,stdin);
- sscanf(line,"%lf",&thresh);
- /* Read and check the parameter values for the tests. */
-
- /* Values of N */
- fgets(line,80,stdin);
- #ifdef USE64BITINT
- sscanf(line,"%ld",&nidim);
- #else
- sscanf(line,"%d",&nidim);
- #endif
-
- if (nidim < 1 || nidim > 9) {
- fprintf(stderr,"NUMBER OF VALUES OF N IS LESS THAN 1 OR GREATER THAN 9");
- goto L220;
- }
- fgets(line,80,stdin);
- #ifdef USE64BITINT
- sscanf(line,"%ld %ld %ld %ld %ld %ld %ld %ld %ld",&idim[0],&idim[1],&idim[2],
- &idim[3],&idim[4],&idim[5],&idim[6],&idim[7],&idim[8]);
- #else
- sscanf(line,"%d %d %d %d %d %d %d %d %d",&idim[0],&idim[1],&idim[2],
- &idim[3],&idim[4],&idim[5],&idim[6],&idim[7],&idim[8]);
- #endif
- i__1 = nidim;
- for (i__ = 1; i__ <= i__1; ++i__) {
- if (idim[i__ - 1] < 0 || idim[i__ - 1] > 65) {
- fprintf(stderr,"VALUE OF N IS LESS THAN 0 OR GREATER THAN 65\n");
- goto L220;
- }
- /* L10: */
- }
- /* Values of ALPHA */
- fgets(line,80,stdin);
- #ifdef USE64BITINT
- sscanf(line,"%ld",&nalf);
- #else
- sscanf(line,"%d",&nalf);
- #endif
- if (nalf < 1 || nalf > 7) {
- fprintf(stderr,"VALUE OF ALPHA IS LESS THAN 0 OR GREATER THAN 7\n");
- goto L220;
- }
- fgets(line,80,stdin);
- sscanf(line,"%lf %lf %lf %lf %lf %lf %lf",&alf[0],&alf[1],&alf[2],&alf[3],&alf[4],&alf[5],&alf[6]);
-
- /* Values of BETA */
- fgets(line,80,stdin);
- #ifdef USE64BITINT
- sscanf(line,"%ld",&nbet);
- #else
- sscanf(line,"%d",&nbet);
- #endif
- if (nalf < 1 || nbet > 7) {
- fprintf(stderr,"VALUE OF BETA IS LESS THAN 0 OR GREATER THAN 7\n");
- goto L220;
- }
- fgets(line,80,stdin);
- sscanf(line,"%lf %lf %lf %lf %lf %lf %lf",&bet[0],&bet[1],&bet[2],&bet[3],&bet[4],&bet[5],&bet[6]);
-
- /* Report values of parameters. */
-
- printf("TESTS OF THE DOUBLE PRECISION LEVEL 3 BLAS\nTHE FOLLOWING PARAMETER VALUES WILL BE USED:\n");
- printf(" FOR N");
- for (i__ =1; i__ <=nidim;++i__) printf(" %d",idim[i__-1]);
- printf("\n");
- printf(" FOR ALPHA");
- for (i__ =1; i__ <=nalf;++i__) printf(" %f",alf[i__-1]);
- printf("\n");
- printf(" FOR BETA");
- for (i__ =1; i__ <=nbet;++i__) printf(" %f",bet[i__-1]);
- printf("\n");
-
- if (! tsterr) {
- printf(" ERROR-EXITS WILL NOT BE TESTED\n");
- }
- printf("ROUTINES PASS COMPUTATIONAL TESTS IF TEST RATIO IS LESS THAN %f\n",thresh);
-
- rorder = FALSE_;
- corder = FALSE_;
- if (layout == 2) {
- rorder = TRUE_;
- corder = TRUE_;
- printf("COLUMN-MAJOR AND ROW-MAJOR DATA LAYOUTS ARE TESTED\n");
- } else if (layout == 1) {
- rorder = TRUE_;
- printf("ROW-MAJOR DATA LAYOUT IS TESTED\n");
- } else if (layout == 0) {
- corder = TRUE_;
- printf("COLUMN-MAJOR DATA LAYOUT IS TESTED\n");
- }
-
- /* Read names of subroutines and flags which indicate */
- /* whether they are to be tested. */
-
- for (i__ = 1; i__ <= 6; ++i__) {
- ltest[i__ - 1] = FALSE_;
- /* L20: */
- }
- L30:
- if (! fgets(line,80,stdin)) {
- goto L60;
- }
- i__1 = sscanf(line,"%12c %c",snamet,&tmpchar);
- ltestt=FALSE_;
- if (tmpchar=='T')ltestt=TRUE_;
- if (i__1 < 2) {
- goto L60;
- }
- for (i__ = 1; i__ <= 6; ++i__) {
- if (s_cmp(snamet, snames[i__ - 1] , (ftnlen)12, (ftnlen)12) ==
- 0) {
- goto L50;
- }
- /* L40: */
- }
- printf("SUBPROGRAM NAME %s NOT RECOGNIZED\n****** TESTS ABANDONED ******\n",snamet);
- exit(1);
-
-
- L50:
- ltest[i__ - 1] = ltestt;
- goto L30;
-
- L60:
- /* cl__1.cerr = 0;
- cl__1.cunit = 5;
- cl__1.csta = 0;
- f_clos(&cl__1);*/
-
- /* Compute EPS (the machine precision). */
-
- eps = 1.;
- L70:
- d__1 = eps + 1.;
- if (ddiff_(&d__1, &c_b90) == 0.) {
- goto L80;
- }
- eps *= .5;
- goto L70;
- L80:
- eps += eps;
- printf("RELATIVE MACHINE PRECISION IS TAKEN TO BE %9.1g\n",eps);
-
- /* Check the reliability of DMMCH using exact data. */
-
- n = 32;
- i__1 = n;
- for (j = 1; j <= i__1; ++j) {
- i__2 = n;
- for (i__ = 1; i__ <= i__2; ++i__) {
- /* Computing MAX */
- i__3 = i__ - j + 1;
- ab[i__ + j * 65 - 66] = (doublereal) f2cmax(i__3,0);
- /* L90: */
- }
- ab[j + 4224] = (doublereal) j;
- ab[(j + 65) * 65 - 65] = (doublereal) j;
- c__[j - 1] = 0.;
- /* L100: */
- }
- i__1 = n;
- for (j = 1; j <= i__1; ++j) {
- cc[j - 1] = (doublereal) (j * ((j + 1) * j) / 2 - (j + 1) * j * (j -
- 1) / 3);
- /* L110: */
- }
- /* CC holds the exact result. On exit from DMMCH CT holds */
- /* the result computed by DMMCH. */
- *(unsigned char *)transa = 'N';
- *(unsigned char *)transb = 'N';
- dmmch_(transa, transb, &n, &c__1, &n, &c_b90, ab, &c__65, &ab[4225], &
- c__65, &c_b104, c__, &c__65, ct, g, cc, &c__65, &eps, &err, &
- fatal, &c__6, &c_true, (ftnlen)1, (ftnlen)1);
- same = lde_(cc, ct, &n);
- if (! same || err != 0.) {
- printf("ERROR IN DMMCH - IN-LINE DOT PRODUCTS ARE BEING EVALUATED WRONGLY\n");
- printf("DMMCH WAS CALLED WITH TRANSA = %s AND TRANSB = %s\n", transa,transb);
- printf("AND RETURNED SAME = %c AND ERR = %12.3f.\n",(same==FALSE_? 'F':'T'),err);
- printf("THIS MAY BE DUE TO FAULTS IN THE ARITHMETIC OR THE COMPILER.\n");
- printf("****** TESTS ABANDONED ******\n");
- exit(1);
- }
- *(unsigned char *)transb = 'T';
- dmmch_(transa, transb, &n, &c__1, &n, &c_b90, ab, &c__65, &ab[4225], &
- c__65, &c_b104, c__, &c__65, ct, g, cc, &c__65, &eps, &err, &
- fatal, &c__6, &c_true, (ftnlen)1, (ftnlen)1);
- same = lde_(cc, ct, &n);
- if (! same || err != 0.) {
- printf("ERROR IN DMMCH - IN-LINE DOT PRODUCTS ARE BEING EVALUATED WRONGLY\n");
- printf("DMMCH WAS CALLED WITH TRANSA = %s AND TRANSB = %s\n", transa,transb);
- printf("AND RETURNED SAME = %c AND ERR = %12.3f.\n",(same==FALSE_? 'F':'T'),err);
- printf("THIS MAY BE DUE TO FAULTS IN THE ARITHMETIC OR THE COMPILER.\n");
- printf("****** TESTS ABANDONED ******\n");
- exit(1);
- }
- i__1 = n;
- for (j = 1; j <= i__1; ++j) {
- ab[j + 4224] = (doublereal) (n - j + 1);
- ab[(j + 65) * 65 - 65] = (doublereal) (n - j + 1);
- /* L120: */
- }
- i__1 = n;
- for (j = 1; j <= i__1; ++j) {
- cc[n - j] = (doublereal) (j * ((j + 1) * j) / 2 - (j + 1) * j * (j -
- 1) / 3);
- /* L130: */
- }
- *(unsigned char *)transa = 'T';
- *(unsigned char *)transb = 'N';
- dmmch_(transa, transb, &n, &c__1, &n, &c_b90, ab, &c__65, &ab[4225], &
- c__65, &c_b104, c__, &c__65, ct, g, cc, &c__65, &eps, &err, &
- fatal, &c__6, &c_true, (ftnlen)1, (ftnlen)1);
- same = lde_(cc, ct, &n);
- if (! same || err != 0.) {
- printf("ERROR IN DMMCH - IN-LINE DOT PRODUCTS ARE BEING EVALUATED WRONGLY\n");
- printf("DMMCH WAS CALLED WITH TRANSA = %s AND TRANSB = %s\n", transa,transb);
- printf("AND RETURNED SAME = %c AND ERR = %12.3f.\n",(same==FALSE_? 'F':'T'),err);
- printf("THIS MAY BE DUE TO FAULTS IN THE ARITHMETIC OR THE COMPILER.\n");
- printf("****** TESTS ABANDONED ******\n");
- exit(1);
- }
- *(unsigned char *)transb = 'T';
- dmmch_(transa, transb, &n, &c__1, &n, &c_b90, ab, &c__65, &ab[4225], &
- c__65, &c_b104, c__, &c__65, ct, g, cc, &c__65, &eps, &err, &
- fatal, &c__6, &c_true, (ftnlen)1, (ftnlen)1);
- same = lde_(cc, ct, &n);
- if (! same || err != 0.) {
- }
-
- /* Test each subroutine in turn. */
-
- for (isnum = 1; isnum <= 6; ++isnum) {
- if (! ltest[isnum - 1]) {
- /* Subprogram is not to be tested. */
- printf("%12s WAS NOT TESTED\n",snames[isnum-1]);
- } else {
- s_copy(srnamc_1.srnamt, snames[isnum - 1], (ftnlen)12, (
- ftnlen)12);
- /* Test error exits. */
- if (tsterr) {
- cd3chke_(snames[isnum - 1], (ftnlen)12);
- }
- /* Test computations. */
- infoc_1.infot = 0;
- infoc_1.ok = TRUE_;
- fatal = FALSE_;
- switch ((int)isnum) {
- case 1: goto L140;
- case 2: goto L150;
- case 3: goto L160;
- case 4: goto L160;
- case 5: goto L170;
- case 6: goto L180;
- }
- /* Test DGEMM, 01. */
- L140:
- if (corder) {
- dchk1_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, &ab[4225], bb, bs, c__,
- cc, cs, ct, g, &c__0, (ftnlen)12);
- }
- if (rorder) {
- dchk1_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, &ab[4225], bb, bs, c__,
- cc, cs, ct, g, &c__1, (ftnlen)12);
- }
- goto L190;
- /* Test DSYMM, 02. */
- L150:
- if (corder) {
- dchk2_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, &ab[4225], bb, bs, c__,
- cc, cs, ct, g, &c__0, (ftnlen)12);
- }
- if (rorder) {
- dchk2_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, &ab[4225], bb, bs, c__,
- cc, cs, ct, g, &c__1, (ftnlen)12);
- }
- goto L190;
- /* Test DTRMM, 03, DTRSM, 04. */
- L160:
- if (corder) {
- dchk3_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- c__65, ab, aa, as, &ab[4225], bb, bs, ct, g, c__, &
- c__0, (ftnlen)12);
- }
- if (rorder) {
- dchk3_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- c__65, ab, aa, as, &ab[4225], bb, bs, ct, g, c__, &
- c__1, (ftnlen)12);
- }
- goto L190;
- /* Test DSYRK, 05. */
- L170:
- if (corder) {
- dchk4_(snames[isnum -1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, &ab[4225], bb, bs, c__,
- cc, cs, ct, g, &c__0, (ftnlen)12);
- }
- if (rorder) {
- dchk4_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, &ab[4225], bb, bs, c__,
- cc, cs, ct, g, &c__1, (ftnlen)12);
- }
- goto L190;
- /* Test DSYR2K, 06. */
- L180:
- if (corder) {
- dchk5_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, bb, bs, c__, cc, cs,
- ct, g, w, &c__0, (ftnlen)12);
- }
- if (rorder) {
- dchk5_(snames[isnum - 1], &eps, &thresh, &c__6, &ntra,
- &trace, &rewi, &fatal, &nidim, idim, &nalf, alf, &
- nbet, bet, &c__65, ab, aa, as, bb, bs, c__, cc, cs,
- ct, g, w, &c__1, (ftnlen)12);
- }
- goto L190;
-
- L190:
- if (fatal && sfatal) {
- goto L210;
- }
- }
- /* L200: */
- }
- printf("\nEND OF TESTS\n");
- goto L230;
-
- L210:
- printf("\n****** FATAL ERROR - TESTS ABANDONED ******\n");
- goto L230;
-
- L220:
- printf("AMEND DATA FILE OR INCREASE ARRAY SIZES IN PROGRAM\n");
- printf("****** TESTS ABANDONED ******\n");
-
- L230:
- if (trace) {
- /* cl__1.cerr = 0;
- cl__1.cunit = ntra;
- cl__1.csta = 0;
- f_clos(&cl__1);*/
- }
- /* cl__1.cerr = 0;
- cl__1.cunit = 6;
- cl__1.csta = 0;
- f_clos(&cl__1);*/
- exit(0);
-
- /* End of DBLAT3. */
-
- } /* MAIN__ */
-
- /* Subroutine */ int dchk1_(sname, eps, thresh, nout, ntra, trace, rewi,
- fatal, nidim, idim, nalf, alf, nbet, bet, nmax, a, aa, as, b, bb, bs,
- c__, cc, cs, ct, g, iorder, sname_len)
- char *sname;
- doublereal *eps, *thresh;
- integer *nout, *ntra;
- logical *trace, *rewi, *fatal;
- integer *nidim, *idim, *nalf;
- doublereal *alf;
- integer *nbet;
- doublereal *bet;
- integer *nmax;
- doublereal *a, *aa, *as, *b, *bb, *bs, *c__, *cc, *cs, *ct, *g;
- integer *iorder;
- ftnlen sname_len;
- {
- /* Initialized data */
-
- static char ich[3+1] = "NTC";
-
- /* System generated locals */
- integer a_dim1, a_offset, b_dim1, b_offset, c_dim1, c_offset, i__1, i__2,
- i__3, i__4, i__5, i__6;
-
- /* Builtin functions */
- integer f_rew(), s_wsfe(), e_wsfe(), do_fio();
-
- /* Local variables */
- static doublereal beta;
- static integer ldas, ldbs, ldcs;
- static logical same, null;
- static integer i__, k, m, n;
- extern /* Subroutine */ int dmake_();
- static doublereal alpha;
- extern /* Subroutine */ int dmmch_();
- static logical isame[13], trana, tranb;
- static integer nargs;
- static logical reset;
- extern /* Subroutine */ void dprcn1_();
- static integer ia, ib, ma, mb, na, nb, nc, ik, im, in;
- extern /* Subroutine */ int cdgemm_();
- static integer ks, ms, ns;
- extern logical lderes_();
- static char tranas[1], tranbs[1], transa[1], transb[1];
- static doublereal errmax;
- static integer ica, icb, laa, lbb, lda, lcc, ldb, ldc;
- extern logical lde_();
- static doublereal als, bls, err;
-
- /* Tests DGEMM. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Local Arrays .. */
- /* .. External Functions .. */
- /* .. External Subroutines .. */
- /* .. Intrinsic Functions .. */
- /* .. Scalars in Common .. */
- /* .. Common blocks .. */
- /* .. Data statements .. */
- /* Parameter adjustments */
- --idim;
- --alf;
- --bet;
- --g;
- --ct;
- --cs;
- --cc;
- c_dim1 = *nmax;
- c_offset = 1 + c_dim1 * 1;
- c__ -= c_offset;
- --bs;
- --bb;
- b_dim1 = *nmax;
- b_offset = 1 + b_dim1 * 1;
- b -= b_offset;
- --as;
- --aa;
- a_dim1 = *nmax;
- a_offset = 1 + a_dim1 * 1;
- a -= a_offset;
-
- /* Function Body */
- /* .. Executable Statements .. */
-
- nargs = 13;
- nc = 0;
- reset = TRUE_;
- errmax = 0.;
-
- i__1 = *nidim;
- for (im = 1; im <= i__1; ++im) {
- m = idim[im];
-
- i__2 = *nidim;
- for (in = 1; in <= i__2; ++in) {
- n = idim[in];
- /* Set LDC to 1 more than minimum value if room. */
- ldc = m;
- if (ldc < *nmax) {
- ++ldc;
- }
- /* Skip tests if not enough room. */
- if (ldc > *nmax) {
- goto L100;
- }
- lcc = ldc * n;
- null = n <= 0 || m <= 0;
-
- i__3 = *nidim;
- for (ik = 1; ik <= i__3; ++ik) {
- k = idim[ik];
-
- for (ica = 1; ica <= 3; ++ica) {
- *(unsigned char *)transa = *(unsigned char *)&ich[ica - 1]
- ;
- trana = *(unsigned char *)transa == 'T' || *(unsigned
- char *)transa == 'C';
-
- if (trana) {
- ma = k;
- na = m;
- } else {
- ma = m;
- na = k;
- }
- /* Set LDA to 1 more than minimum value if room. */
- lda = ma;
- if (lda < *nmax) {
- ++lda;
- }
- /* Skip tests if not enough room. */
- if (lda > *nmax) {
- goto L80;
- }
- laa = lda * na;
-
- /* Generate the matrix A. */
-
- dmake_("GE", " ", " ", &ma, &na, &a[a_offset], nmax, &aa[
- 1], &lda, &reset, &c_b104, (ftnlen)2, (ftnlen)1, (
- ftnlen)1);
-
- for (icb = 1; icb <= 3; ++icb) {
- *(unsigned char *)transb = *(unsigned char *)&ich[icb
- - 1];
- tranb = *(unsigned char *)transb == 'T' || *(unsigned
- char *)transb == 'C';
-
- if (tranb) {
- mb = n;
- nb = k;
- } else {
- mb = k;
- nb = n;
- }
- /* Set LDB to 1 more than minimum value if room. */
- ldb = mb;
- if (ldb < *nmax) {
- ++ldb;
- }
- /* Skip tests if not enough room. */
- if (ldb > *nmax) {
- goto L70;
- }
- lbb = ldb * nb;
-
- /* Generate the matrix B. */
-
- dmake_("GE", " ", " ", &mb, &nb, &b[b_offset], nmax, &
- bb[1], &ldb, &reset, &c_b104, (ftnlen)2, (
- ftnlen)1, (ftnlen)1);
-
- i__4 = *nalf;
- for (ia = 1; ia <= i__4; ++ia) {
- alpha = alf[ia];
-
- i__5 = *nbet;
- for (ib = 1; ib <= i__5; ++ib) {
- beta = bet[ib];
-
- /* Generate the matrix C. */
-
- dmake_("GE", " ", " ", &m, &n, &c__[c_offset],
- nmax, &cc[1], &ldc, &reset, &c_b104,
- (ftnlen)2, (ftnlen)1, (ftnlen)1);
-
- ++nc;
-
- /* Save every datum before calling the */
- /* subroutine. */
-
- *(unsigned char *)tranas = *(unsigned char *)
- transa;
- *(unsigned char *)tranbs = *(unsigned char *)
- transb;
- ms = m;
- ns = n;
- ks = k;
- als = alpha;
- i__6 = laa;
- for (i__ = 1; i__ <= i__6; ++i__) {
- as[i__] = aa[i__];
- /* L10: */
- }
- ldas = lda;
- i__6 = lbb;
- for (i__ = 1; i__ <= i__6; ++i__) {
- bs[i__] = bb[i__];
- /* L20: */
- }
- ldbs = ldb;
- bls = beta;
- i__6 = lcc;
- for (i__ = 1; i__ <= i__6; ++i__) {
- cs[i__] = cc[i__];
- /* L30: */
- }
- ldcs = ldc;
-
- /* Call the subroutine. */
-
- if (*trace) {
- dprcn1_(ntra, &nc, sname, iorder, transa,
- transb, &m, &n, &k, &alpha, &lda,
- &ldb, &beta, &ldc, (ftnlen)12, (
- ftnlen)1, (ftnlen)1);
- }
- if (*rewi) {
- /* al__1.aerr = 0;
- al__1.aunit = *ntra;
- f_rew(&al__1);*/
- }
- cdgemm_(iorder, transa, transb, &m, &n, &k, &
- alpha, &aa[1], &lda, &bb[1], &ldb, &
- beta, &cc[1], &ldc, (ftnlen)1, (
- ftnlen)1);
-
- /* Check if error-exit was taken incorrectly. */
-
- if (! infoc_1.ok) {
- printf(" *** FATAL ERROR - ERROR-CALL MYEXIT TAKEN ON VALID CALL\n");
- *fatal = TRUE_;
- goto L120;
- }
-
- /* See what data changed inside subroutines. */
-
- isame[0] = *(unsigned char *)transa == *(
- unsigned char *)tranas;
- isame[1] = *(unsigned char *)transb == *(
- unsigned char *)tranbs;
- isame[2] = ms == m;
- isame[3] = ns == n;
- isame[4] = ks == k;
- isame[5] = als == alpha;
- isame[6] = lde_(&as[1], &aa[1], &laa);
- isame[7] = ldas == lda;
- isame[8] = lde_(&bs[1], &bb[1], &lbb);
- isame[9] = ldbs == ldb;
- isame[10] = bls == beta;
- if (null) {
- isame[11] = lde_(&cs[1], &cc[1], &lcc);
- } else {
- isame[11] = lderes_("GE", " ", &m, &n, &
- cs[1], &cc[1], &ldc, (ftnlen)2, (
- ftnlen)1);
- }
- isame[12] = ldcs == ldc;
-
- /* If data was incorrectly changed, report */
- /* and return. */
-
- same = TRUE_;
- i__6 = nargs;
- for (i__ = 1; i__ <= i__6; ++i__) {
- same = same && isame[i__ - 1];
- if (! isame[i__ - 1]) {
- printf(" ******* FATAL ERROR - PARAMETER NUMBER %d WAS CHANGED INCORRECTLY *******\n",i__);
- }
- /* L40: */
- }
- if (! same) {
- *fatal = TRUE_;
- goto L120;
- }
-
- if (! null) {
-
- /* Check the result. */
-
- dmmch_(transa, transb, &m, &n, &k, &alpha,
- &a[a_offset], nmax, &b[b_offset],
- nmax, &beta, &c__[c_offset],
- nmax, &ct[1], &g[1], &cc[1], &ldc,
- eps, &err, fatal, nout, &c_true,
- (ftnlen)1, (ftnlen)1);
- errmax = f2cmax(errmax,err);
- /* If got really bad answer, report and */
- /* return. */
- if (*fatal) {
- goto L120;
- }
- }
-
- /* L50: */
- }
-
- /* L60: */
- }
-
- L70:
- ;
- }
-
- L80:
- ;
- }
-
- /* L90: */
- }
-
- L100:
- ;
- }
-
- /* L110: */
- }
-
- /* Report result. */
-
- if (errmax < *thresh) {
- if (*iorder == 0) {
- printf("%s PASSED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- if (*iorder == 1) {
- printf("%s PASSED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- } else {
- if (*iorder == 0) {
- printf("%s COMPLETED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- if (*iorder == 1) {
- printf("%s COMPLETED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- }
- goto L130;
-
- L120:
- printf(" ******* %s FAILED ON CALL NUMBER:\n",sname);
- dprcn1_(nout, &nc, sname, iorder, transa, transb, &m, &n, &k, &alpha, &
- lda, &ldb, &beta, &ldc, (ftnlen)12, (ftnlen)1, (ftnlen)1);
-
- L130:
- return 0;
-
- /* 9995 FORMAT( 1X, I6, ': ', A12,'(''', A1, ''',''', A1, ''',', */
- /* $ 3( I3, ',' ), F4.1, ', A,', I3, ', B,', I3, ',', F4.1, ', ', */
- /* $ 'C,', I3, ').' ) */
-
- /* End of DCHK1. */
-
- } /* dchk1_ */
-
- /* Subroutine */ void dprcn1_(nout, nc, sname, iorder, transa, transb, m, n, k,
- alpha, lda, ldb, beta, ldc, sname_len, transa_len, transb_len)
- integer *nout, *nc;
- char *sname;
- integer *iorder;
- char *transa, *transb;
- integer *m, *n, *k;
- doublereal *alpha;
- integer *lda, *ldb;
- doublereal *beta;
- integer *ldc;
- ftnlen sname_len;
- ftnlen transa_len;
- ftnlen transb_len;
- {
- /* Builtin functions */
- integer s_wsfe(), do_fio(), e_wsfe();
-
- /* Local variables */
- static char crc[14], cta[14], ctb[14];
-
- if (*(unsigned char *)transa == 'N') {
- s_copy(cta, " CblasNoTrans", (ftnlen)14, (ftnlen)14);
- } else if (*(unsigned char *)transa == 'T') {
- s_copy(cta, " CblasTrans", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cta, "CblasConjTrans", (ftnlen)14, (ftnlen)14);
- }
- if (*(unsigned char *)transb == 'N') {
- s_copy(ctb, " CblasNoTrans", (ftnlen)14, (ftnlen)14);
- } else if (*(unsigned char *)transb == 'T') {
- s_copy(ctb, " CblasTrans", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(ctb, "CblasConjTrans", (ftnlen)14, (ftnlen)14);
- }
- if (*iorder == 1) {
- s_copy(crc, " CblasRowMajor", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(crc, " CblasColMajor", (ftnlen)14, (ftnlen)14);
- }
- printf("%6d: %s %s %s %s\n",*nc,sname,crc,cta,ctb);
- printf("%d %d %d %4.1f A, %d, B, %d, %4.1f, C, %d.\n",*m,*n,*k,*alpha,*lda,*ldb,*beta,*ldc);
- } /* dprcn1_ */
-
-
- /* Subroutine */ int dchk2_(sname, eps, thresh, nout, ntra, trace, rewi,
- fatal, nidim, idim, nalf, alf, nbet, bet, nmax, a, aa, as, b, bb, bs,
- c__, cc, cs, ct, g, iorder, sname_len)
- char *sname;
- doublereal *eps, *thresh;
- integer *nout, *ntra;
- logical *trace, *rewi, *fatal;
- integer *nidim, *idim, *nalf;
- doublereal *alf;
- integer *nbet;
- doublereal *bet;
- integer *nmax;
- doublereal *a, *aa, *as, *b, *bb, *bs, *c__, *cc, *cs, *ct, *g;
- integer *iorder;
- ftnlen sname_len;
- {
- /* Initialized data */
-
- static char ichs[2+1] = "LR";
- static char ichu[2+1] = "UL";
-
- /* System generated locals */
- integer a_dim1, a_offset, b_dim1, b_offset, c_dim1, c_offset, i__1, i__2,
- i__3, i__4, i__5;
-
- /* Builtin functions */
- integer f_rew(), s_wsfe(), e_wsfe(), do_fio();
-
- /* Local variables */
- static doublereal beta;
- static integer ldas, ldbs, ldcs;
- static logical same;
- static char side[1];
- static logical left, null;
- static char uplo[1];
- static integer i__, m, n;
- extern /* Subroutine */ int dmake_();
- static doublereal alpha;
- extern /* Subroutine */ int dmmch_();
- static logical isame[13];
- static char sides[1];
- static integer nargs;
- static logical reset;
- static char uplos[1];
- extern /* Subroutine */ void dprcn2_();
- static integer ia, ib, na, nc, im, in, ms, ns;
- extern logical lderes_();
- extern /* Subroutine */ int cdsymm_();
- static doublereal errmax;
- static integer laa, lbb, lda, lcc, ldb, ldc;
- extern logical lde_();
- static integer ics;
- static doublereal als, bls;
- static integer icu;
- static doublereal err;
-
- /* Tests DSYMM. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Local Arrays .. */
- /* .. External Functions .. */
- /* .. External Subroutines .. */
- /* .. Intrinsic Functions .. */
- /* .. Scalars in Common .. */
- /* .. Common blocks .. */
- /* .. Data statements .. */
- /* Parameter adjustments */
- --idim;
- --alf;
- --bet;
- --g;
- --ct;
- --cs;
- --cc;
- c_dim1 = *nmax;
- c_offset = 1 + c_dim1 * 1;
- c__ -= c_offset;
- --bs;
- --bb;
- b_dim1 = *nmax;
- b_offset = 1 + b_dim1 * 1;
- b -= b_offset;
- --as;
- --aa;
- a_dim1 = *nmax;
- a_offset = 1 + a_dim1 * 1;
- a -= a_offset;
-
- /* Function Body */
- /* .. Executable Statements .. */
-
- nargs = 12;
- nc = 0;
- reset = TRUE_;
- errmax = 0.;
-
- i__1 = *nidim;
- for (im = 1; im <= i__1; ++im) {
- m = idim[im];
-
- i__2 = *nidim;
- for (in = 1; in <= i__2; ++in) {
- n = idim[in];
- /* Set LDC to 1 more than minimum value if room. */
- ldc = m;
- if (ldc < *nmax) {
- ++ldc;
- }
- /* Skip tests if not enough room. */
- if (ldc > *nmax) {
- goto L90;
- }
- lcc = ldc * n;
- null = n <= 0 || m <= 0;
-
- /* Set LDB to 1 more than minimum value if room. */
- ldb = m;
- if (ldb < *nmax) {
- ++ldb;
- }
- /* Skip tests if not enough room. */
- if (ldb > *nmax) {
- goto L90;
- }
- lbb = ldb * n;
-
- /* Generate the matrix B. */
-
- dmake_("GE", " ", " ", &m, &n, &b[b_offset], nmax, &bb[1], &ldb, &
- reset, &c_b104, (ftnlen)2, (ftnlen)1, (ftnlen)1);
-
- for (ics = 1; ics <= 2; ++ics) {
- *(unsigned char *)side = *(unsigned char *)&ichs[ics - 1];
- left = *(unsigned char *)side == 'L';
-
- if (left) {
- na = m;
- } else {
- na = n;
- }
- /* Set LDA to 1 more than minimum value if room. */
- lda = na;
- if (lda < *nmax) {
- ++lda;
- }
- /* Skip tests if not enough room. */
- if (lda > *nmax) {
- goto L80;
- }
- laa = lda * na;
-
- for (icu = 1; icu <= 2; ++icu) {
- *(unsigned char *)uplo = *(unsigned char *)&ichu[icu - 1];
-
- /* Generate the symmetric matrix A. */
-
- dmake_("SY", uplo, " ", &na, &na, &a[a_offset], nmax, &aa[
- 1], &lda, &reset, &c_b104, (ftnlen)2, (ftnlen)1, (
- ftnlen)1);
-
- i__3 = *nalf;
- for (ia = 1; ia <= i__3; ++ia) {
- alpha = alf[ia];
-
- i__4 = *nbet;
- for (ib = 1; ib <= i__4; ++ib) {
- beta = bet[ib];
-
- /* Generate the matrix C. */
-
- dmake_("GE", " ", " ", &m, &n, &c__[c_offset],
- nmax, &cc[1], &ldc, &reset, &c_b104, (
- ftnlen)2, (ftnlen)1, (ftnlen)1);
-
- ++nc;
-
- /* Save every datum before calling the */
- /* subroutine. */
-
- *(unsigned char *)sides = *(unsigned char *)side;
- *(unsigned char *)uplos = *(unsigned char *)uplo;
- ms = m;
- ns = n;
- als = alpha;
- i__5 = laa;
- for (i__ = 1; i__ <= i__5; ++i__) {
- as[i__] = aa[i__];
- /* L10: */
- }
- ldas = lda;
- i__5 = lbb;
- for (i__ = 1; i__ <= i__5; ++i__) {
- bs[i__] = bb[i__];
- /* L20: */
- }
- ldbs = ldb;
- bls = beta;
- i__5 = lcc;
- for (i__ = 1; i__ <= i__5; ++i__) {
- cs[i__] = cc[i__];
- /* L30: */
- }
- ldcs = ldc;
-
- /* Call the subroutine. */
-
- if (*trace) {
- dprcn2_(ntra, &nc, sname, iorder, side, uplo,
- &m, &n, &alpha, &lda, &ldb, &beta, &
- ldc, (ftnlen)12, (ftnlen)1, (ftnlen)1)
- ;
- }
- if (*rewi) {
- /* al__1.aerr = 0;
- al__1.aunit = *ntra;
- f_rew(&al__1);*/
- }
- cdsymm_(iorder, side, uplo, &m, &n, &alpha, &aa[1]
- , &lda, &bb[1], &ldb, &beta, &cc[1], &ldc,
- (ftnlen)1, (ftnlen)1);
-
- /* Check if error-exit was taken incorrectly. */
-
- if (! infoc_1.ok) {
- printf("*** FATAL ERROR - ERROR-CALL MYEXIT TAKEN ON VALID CALL\n");
- *fatal = TRUE_;
- goto L110;
- }
-
- /* See what data changed inside subroutines. */
-
- isame[0] = *(unsigned char *)sides == *(unsigned
- char *)side;
- isame[1] = *(unsigned char *)uplos == *(unsigned
- char *)uplo;
- isame[2] = ms == m;
- isame[3] = ns == n;
- isame[4] = als == alpha;
- isame[5] = lde_(&as[1], &aa[1], &laa);
- isame[6] = ldas == lda;
- isame[7] = lde_(&bs[1], &bb[1], &lbb);
- isame[8] = ldbs == ldb;
- isame[9] = bls == beta;
- if (null) {
- isame[10] = lde_(&cs[1], &cc[1], &lcc);
- } else {
- isame[10] = lderes_("GE", " ", &m, &n, &cs[1],
- &cc[1], &ldc, (ftnlen)2, (ftnlen)1);
- }
- isame[11] = ldcs == ldc;
-
- /* If data was incorrectly changed, report and */
- /* return. */
-
- same = TRUE_;
- i__5 = nargs;
- for (i__ = 1; i__ <= i__5; ++i__) {
- same = same && isame[i__ - 1];
- if (! isame[i__ - 1]) {
- printf(" ******* FATAL ERROR - PARAMETER NUMBER %d WAS CHANGED INCORRECTLY *******\n",i__);
- }
- /* L40: */
- }
- if (! same) {
- *fatal = TRUE_;
- goto L110;
- }
-
- if (! null) {
-
- /* Check the result. */
-
- if (left) {
- dmmch_("N", "N", &m, &n, &m, &alpha, &a[
- a_offset], nmax, &b[b_offset],
- nmax, &beta, &c__[c_offset], nmax,
- &ct[1], &g[1], &cc[1], &ldc, eps,
- &err, fatal, nout, &c_true, (
- ftnlen)1, (ftnlen)1);
- } else {
- dmmch_("N", "N", &m, &n, &n, &alpha, &b[
- b_offset], nmax, &a[a_offset],
- nmax, &beta, &c__[c_offset], nmax,
- &ct[1], &g[1], &cc[1], &ldc, eps,
- &err, fatal, nout, &c_true, (
- ftnlen)1, (ftnlen)1);
- }
- errmax = f2cmax(errmax,err);
- /* If got really bad answer, report and */
- /* return. */
- if (*fatal) {
- goto L110;
- }
- }
-
- /* L50: */
- }
-
- /* L60: */
- }
-
- /* L70: */
- }
-
- L80:
- ;
- }
-
- L90:
- ;
- }
-
- /* L100: */
- }
-
- /* Report result. */
-
- if (errmax < *thresh) {
- if (*iorder == 0) {
- printf("%s PASSED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- if (*iorder == 1) {
- printf("%s PASSED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- } else {
- if (*iorder == 0) {
- printf("%s COMPLETED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- if (*iorder == 1) {
- printf("%s COMPLETED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- }
- goto L120;
-
- L110:
- printf(" ******* %s FAILED ON CALL NUMBER:\n",sname);
- dprcn2_(nout, &nc, sname, iorder, side, uplo, &m, &n, &alpha, &lda, &ldb,
- &beta, &ldc, (ftnlen)12, (ftnlen)1, (ftnlen)1);
-
- L120:
- return 0;
-
- /* 9995 FORMAT( 1X, I6, ': ', A12,'(', 2( '''', A1, ''',' ), 2( I3, ',' ), */
- /* $ F4.1, ', A,', I3, ', B,', I3, ',', F4.1, ', C,', I3, ') ', */
- /* $ ' .' ) */
-
- /* End of DCHK2. */
-
- } /* dchk2_ */
-
-
- /* Subroutine */ void dprcn2_(nout, nc, sname, iorder, side, uplo, m, n, alpha,
- lda, ldb, beta, ldc, sname_len, side_len, uplo_len)
- integer *nout, *nc;
- char *sname;
- integer *iorder;
- char *side, *uplo;
- integer *m, *n;
- doublereal *alpha;
- integer *lda, *ldb;
- doublereal *beta;
- integer *ldc;
- ftnlen sname_len;
- ftnlen side_len;
- ftnlen uplo_len;
- {
- /* Builtin functions */
- integer s_wsfe(), do_fio(), e_wsfe();
-
- /* Local variables */
- static char cs[14], cu[14], crc[14];
-
- if (*(unsigned char *)side == 'L') {
- s_copy(cs, " CblasLeft", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cs, " CblasRight", (ftnlen)14, (ftnlen)14);
- }
- if (*(unsigned char *)uplo == 'U') {
- s_copy(cu, " CblasUpper", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cu, " CblasLower", (ftnlen)14, (ftnlen)14);
- }
- if (*iorder == 1) {
- s_copy(crc, " CblasRowMajor", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(crc, " CblasColMajor", (ftnlen)14, (ftnlen)14);
- }
- printf("%6d: %s %s %s %s\n",*nc,sname,crc,cs,cu);
- printf("%d %d %4.1f A, %d, B, %d, %4.1f C, %d.\n",*m,*n,*alpha,*lda,*ldb,*beta,*ldc);
- } /* dprcn2_ */
-
-
- /* Subroutine */ int dchk3_(sname, eps, thresh, nout, ntra, trace, rewi,
- fatal, nidim, idim, nalf, alf, nmax, a, aa, as, b, bb, bs, ct, g, c__,
- iorder, sname_len)
- char *sname;
- doublereal *eps, *thresh;
- integer *nout, *ntra;
- logical *trace, *rewi, *fatal;
- integer *nidim, *idim, *nalf;
- doublereal *alf;
- integer *nmax;
- doublereal *a, *aa, *as, *b, *bb, *bs, *ct, *g, *c__;
- integer *iorder;
- ftnlen sname_len;
- {
- /* Initialized data */
-
- static char ichu[2+1] = "UL";
- static char icht[3+1] = "NTC";
- static char ichd[2+1] = "UN";
- static char ichs[2+1] = "LR";
-
- /* System generated locals */
- integer a_dim1, a_offset, b_dim1, b_offset, c_dim1, c_offset, i__1, i__2,
- i__3, i__4, i__5;
-
- /* Local variables */
- static char diag[1];
- static integer ldas, ldbs;
- static logical same;
- static char side[1];
- static logical left, null;
- static char uplo[1];
- static integer i__, j, m, n;
- extern /* Subroutine */ int dmake_();
- static doublereal alpha;
- static char diags[1];
- extern /* Subroutine */ int dmmch_();
- static logical isame[13];
- static char sides[1];
- static integer nargs;
- static logical reset;
- static char uplos[1];
- extern /* Subroutine */ void dprcn3_();
- static integer ia, na, nc, im, in, ms, ns;
- extern logical lderes_();
- extern /* Subroutine */ int cdtrmm_();
- static char tranas[1], transa[1];
- extern /* Subroutine */ int cdtrsm_();
- static doublereal errmax;
- static integer laa, icd, lbb, lda, ldb;
- extern logical lde_();
- static integer ics;
- static doublereal als;
- static integer ict, icu;
- static doublereal err;
-
- /* Tests DTRMM and DTRSM. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Local Arrays .. */
- /* .. External Functions .. */
- /* .. External Subroutines .. */
- /* .. Intrinsic Functions .. */
- /* .. Scalars in Common .. */
- /* .. Common blocks .. */
- /* .. Data statements .. */
- /* Parameter adjustments */
- --idim;
- --alf;
- c_dim1 = *nmax;
- c_offset = 1 + c_dim1 * 1;
- c__ -= c_offset;
- --g;
- --ct;
- --bs;
- --bb;
- b_dim1 = *nmax;
- b_offset = 1 + b_dim1 * 1;
- b -= b_offset;
- --as;
- --aa;
- a_dim1 = *nmax;
- a_offset = 1 + a_dim1 * 1;
- a -= a_offset;
-
- /* Function Body */
- /* .. Executable Statements .. */
-
- nargs = 11;
- nc = 0;
- reset = TRUE_;
- errmax = 0.;
- /* Set up zero matrix for DMMCH. */
- i__1 = *nmax;
- for (j = 1; j <= i__1; ++j) {
- i__2 = *nmax;
- for (i__ = 1; i__ <= i__2; ++i__) {
- c__[i__ + j * c_dim1] = 0.;
- /* L10: */
- }
- /* L20: */
- }
-
- i__1 = *nidim;
- for (im = 1; im <= i__1; ++im) {
- m = idim[im];
-
- i__2 = *nidim;
- for (in = 1; in <= i__2; ++in) {
- n = idim[in];
- /* Set LDB to 1 more than minimum value if room. */
- ldb = m;
- if (ldb < *nmax) {
- ++ldb;
- }
- /* Skip tests if not enough room. */
- if (ldb > *nmax) {
- goto L130;
- }
- lbb = ldb * n;
- null = m <= 0 || n <= 0;
-
- for (ics = 1; ics <= 2; ++ics) {
- *(unsigned char *)side = *(unsigned char *)&ichs[ics - 1];
- left = *(unsigned char *)side == 'L';
- if (left) {
- na = m;
- } else {
- na = n;
- }
- /* Set LDA to 1 more than minimum value if room. */
- lda = na;
- if (lda < *nmax) {
- ++lda;
- }
- /* Skip tests if not enough room. */
- if (lda > *nmax) {
- goto L130;
- }
- laa = lda * na;
-
- for (icu = 1; icu <= 2; ++icu) {
- *(unsigned char *)uplo = *(unsigned char *)&ichu[icu - 1];
-
- for (ict = 1; ict <= 3; ++ict) {
- *(unsigned char *)transa = *(unsigned char *)&icht[
- ict - 1];
-
- for (icd = 1; icd <= 2; ++icd) {
- *(unsigned char *)diag = *(unsigned char *)&ichd[
- icd - 1];
-
- i__3 = *nalf;
- for (ia = 1; ia <= i__3; ++ia) {
- alpha = alf[ia];
-
- /* Generate the matrix A. */
-
- dmake_("TR", uplo, diag, &na, &na, &a[
- a_offset], nmax, &aa[1], &lda, &reset,
- &c_b104, (ftnlen)2, (ftnlen)1, (
- ftnlen)1);
-
- /* Generate the matrix B. */
-
- dmake_("GE", " ", " ", &m, &n, &b[b_offset],
- nmax, &bb[1], &ldb, &reset, &c_b104, (
- ftnlen)2, (ftnlen)1, (ftnlen)1);
-
- ++nc;
-
- /* Save every datum before calling the */
- /* subroutine. */
-
- *(unsigned char *)sides = *(unsigned char *)
- side;
- *(unsigned char *)uplos = *(unsigned char *)
- uplo;
- *(unsigned char *)tranas = *(unsigned char *)
- transa;
- *(unsigned char *)diags = *(unsigned char *)
- diag;
- ms = m;
- ns = n;
- als = alpha;
- i__4 = laa;
- for (i__ = 1; i__ <= i__4; ++i__) {
- as[i__] = aa[i__];
- /* L30: */
- }
- ldas = lda;
- i__4 = lbb;
- for (i__ = 1; i__ <= i__4; ++i__) {
- bs[i__] = bb[i__];
- /* L40: */
- }
- ldbs = ldb;
-
- /* Call the subroutine. */
-
- if (s_cmp(sname + 9, "mm", (ftnlen)2, (ftnlen)
- 2) == 0) {
- if (*trace) {
- dprcn3_(ntra, &nc, sname, iorder,
- side, uplo, transa, diag, &m,
- &n, &alpha, &lda, &ldb, (
- ftnlen)12, (ftnlen)1, (ftnlen)
- 1, (ftnlen)1, (ftnlen)1);
- }
- if (*rewi) {
- /* al__1.aerr = 0;
- al__1.aunit = *ntra;
- f_rew(&al__1);*/
- }
- cdtrmm_(iorder, side, uplo, transa, diag,
- &m, &n, &alpha, &aa[1], &lda, &bb[
- 1], &ldb, (ftnlen)1, (ftnlen)1, (
- ftnlen)1, (ftnlen)1);
- } else if (s_cmp(sname + 9, "sm", (ftnlen)2, (
- ftnlen)2) == 0) {
- if (*trace) {
- dprcn3_(ntra, &nc, sname, iorder,
- side, uplo, transa, diag, &m,
- &n, &alpha, &lda, &ldb, (
- ftnlen)12, (ftnlen)1, (ftnlen)
- 1, (ftnlen)1, (ftnlen)1);
- }
- if (*rewi) {
- /* al__1.aerr = 0;
- al__1.aunit = *ntra;
- f_rew(&al__1);*/
- }
- cdtrsm_(iorder, side, uplo, transa, diag,
- &m, &n, &alpha, &aa[1], &lda, &bb[
- 1], &ldb, (ftnlen)1, (ftnlen)1, (
- ftnlen)1, (ftnlen)1);
- }
-
- /* Check if error-exit was taken incorrectly. */
-
- if (! infoc_1.ok) {
- printf("*** FATAL ERROR - ERROR-CALL MYEXIT TAKEN ON VALID CALL\n");
- *fatal = TRUE_;
- goto L150;
- }
-
- /* See what data changed inside subroutines. */
-
- isame[0] = *(unsigned char *)sides == *(
- unsigned char *)side;
- isame[1] = *(unsigned char *)uplos == *(
- unsigned char *)uplo;
- isame[2] = *(unsigned char *)tranas == *(
- unsigned char *)transa;
- isame[3] = *(unsigned char *)diags == *(
- unsigned char *)diag;
- isame[4] = ms == m;
- isame[5] = ns == n;
- isame[6] = als == alpha;
- isame[7] = lde_(&as[1], &aa[1], &laa);
- isame[8] = ldas == lda;
- if (null) {
- isame[9] = lde_(&bs[1], &bb[1], &lbb);
- } else {
- isame[9] = lderes_("GE", " ", &m, &n, &bs[
- 1], &bb[1], &ldb, (ftnlen)2, (
- ftnlen)1);
- }
- isame[10] = ldbs == ldb;
-
- /* If data was incorrectly changed, report and */
- /* return. */
-
- same = TRUE_;
- i__4 = nargs;
- for (i__ = 1; i__ <= i__4; ++i__) {
- same = same && isame[i__ - 1];
- if (! isame[i__ - 1]) {
- printf(" ******* FATAL ERROR - PARAMETER NUMBER %d WAS CHANGED INCORRECTLY *******\n",i__);
- }
- /* L50: */
- }
- if (! same) {
- *fatal = TRUE_;
- goto L150;
- }
-
- if (! null) {
- if (s_cmp(sname + 9, "mm", (ftnlen)2, (
- ftnlen)2) == 0) {
-
- /* Check the result. */
-
- if (left) {
- dmmch_(transa, "N", &m, &n, &m, &
- alpha, &a[a_offset], nmax,
- &b[b_offset], nmax, &
- c_b104, &c__[c_offset],
- nmax, &ct[1], &g[1], &bb[
- 1], &ldb, eps, &err,
- fatal, nout, &c_true, (
- ftnlen)1, (ftnlen)1);
- } else {
- dmmch_("N", transa, &m, &n, &n, &
- alpha, &b[b_offset], nmax,
- &a[a_offset], nmax, &
- c_b104, &c__[c_offset],
- nmax, &ct[1], &g[1], &bb[
- 1], &ldb, eps, &err,
- fatal, nout, &c_true, (
- ftnlen)1, (ftnlen)1);
- }
- } else if (s_cmp(sname + 9, "sm", (ftnlen)
- 2, (ftnlen)2) == 0) {
-
- /* Compute approximation to original */
- /* matrix. */
-
- i__4 = n;
- for (j = 1; j <= i__4; ++j) {
- i__5 = m;
- for (i__ = 1; i__ <= i__5; ++i__)
- {
- c__[i__ + j * c_dim1] = bb[i__ + (j - 1) * ldb];
- bb[i__ + (j - 1) * ldb] = alpha * b[i__ + j *
- b_dim1];
- /* L60: */
- }
- /* L70: */
- }
-
- if (left) {
- dmmch_(transa, "N", &m, &n, &m, &
- c_b90, &a[a_offset], nmax,
- &c__[c_offset], nmax, &
- c_b104, &b[b_offset],
- nmax, &ct[1], &g[1], &bb[
- 1], &ldb, eps, &err,
- fatal, nout, &c_false, (
- ftnlen)1, (ftnlen)1);
- } else {
- dmmch_("N", transa, &m, &n, &n, &
- c_b90, &c__[c_offset],
- nmax, &a[a_offset], nmax,
- &c_b104, &b[b_offset],
- nmax, &ct[1], &g[1], &bb[
- 1], &ldb, eps, &err,
- fatal, nout, &c_false, (
- ftnlen)1, (ftnlen)1);
- }
- }
- errmax = f2cmax(errmax,err);
- /* If got really bad answer, report and */
- /* return. */
- if (*fatal) {
- goto L150;
- }
- }
-
- /* L80: */
- }
-
- /* L90: */
- }
-
- /* L100: */
- }
-
- /* L110: */
- }
-
- /* L120: */
- }
-
- L130:
- ;
- }
-
- /* L140: */
- }
-
- /* Report result. */
-
- if (errmax < *thresh) {
- if (*iorder == 0) {
- printf("%s PASSED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- if (*iorder == 1) {
- printf("%s PASSED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- } else {
- if (*iorder == 0) {
- printf("%s COMPLETED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- if (*iorder == 1) {
- printf("%s COMPLETED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- }
- goto L160;
-
- L150:
- printf(" ******* %s FAILED ON CALL NUMBER:\n",sname);
- if (*trace) {
- dprcn3_(ntra, &nc, sname, iorder, side, uplo, transa, diag, &m, &n, &
- alpha, &lda, &ldb, (ftnlen)12, (ftnlen)1, (ftnlen)1, (ftnlen)
- 1, (ftnlen)1);
- }
-
- L160:
- return 0;
-
- /* 9995 FORMAT( 1X, I6, ': ', A12,'(', 4( '''', A1, ''',' ), 2( I3, ',' ), */
- /* $ F4.1, ', A,', I3, ', B,', I3, ') .' ) */
-
- /* End of DCHK3. */
-
- } /* dchk3_ */
-
-
- /* Subroutine */ void dprcn3_(nout, nc, sname, iorder, side, uplo, transa,
- diag, m, n, alpha, lda, ldb, sname_len, side_len, uplo_len,
- transa_len, diag_len)
- integer *nout, *nc;
- char *sname;
- integer *iorder;
- char *side, *uplo, *transa, *diag;
- integer *m, *n;
- doublereal *alpha;
- integer *lda, *ldb;
- ftnlen sname_len;
- ftnlen side_len;
- ftnlen uplo_len;
- ftnlen transa_len;
- ftnlen diag_len;
- {
- /* Builtin functions */
- integer s_wsfe(), do_fio(), e_wsfe();
-
- /* Local variables */
- static char ca[14], cd[14], cs[14], cu[14], crc[14];
-
- if (*(unsigned char *)side == 'L') {
- s_copy(cs, " CblasLeft", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cs, " CblasRight", (ftnlen)14, (ftnlen)14);
- }
- if (*(unsigned char *)uplo == 'U') {
- s_copy(cu, " CblasUpper", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cu, " CblasLower", (ftnlen)14, (ftnlen)14);
- }
- if (*(unsigned char *)transa == 'N') {
- s_copy(ca, " CblasNoTrans", (ftnlen)14, (ftnlen)14);
- } else if (*(unsigned char *)transa == 'T') {
- s_copy(ca, " CblasTrans", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(ca, "CblasConjTrans", (ftnlen)14, (ftnlen)14);
- }
- if (*(unsigned char *)diag == 'N') {
- s_copy(cd, " CblasNonUnit", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cd, " CblasUnit", (ftnlen)14, (ftnlen)14);
- }
- if (*iorder == 1) {
- s_copy(crc, " CblasRowMajor", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(crc, " CblasColMajor", (ftnlen)14, (ftnlen)14);
- }
- printf("%6d: %s %s %s %s\n",*nc,sname,crc,cs,cu);
- printf(" %s %s %d %d %4.1f A %d B %d\n",ca,cd,*m,*n,*alpha,*lda,*ldb);
- } /* dprcn3_ */
-
-
- /* Subroutine */ int dchk4_(sname, eps, thresh, nout, ntra, trace, rewi,
- fatal, nidim, idim, nalf, alf, nbet, bet, nmax, a, aa, as, b, bb, bs,
- c__, cc, cs, ct, g, iorder, sname_len)
- char *sname;
- doublereal *eps, *thresh;
- integer *nout, *ntra;
- logical *trace, *rewi, *fatal;
- integer *nidim, *idim, *nalf;
- doublereal *alf;
- integer *nbet;
- doublereal *bet;
- integer *nmax;
- doublereal *a, *aa, *as, *b, *bb, *bs, *c__, *cc, *cs, *ct, *g;
- integer *iorder;
- ftnlen sname_len;
- {
- /* Initialized data */
-
- static char icht[3+1] = "NTC";
- static char ichu[2+1] = "UL";
-
- /* System generated locals */
- integer a_dim1, a_offset, b_dim1, b_offset, c_dim1, c_offset, i__1, i__2,
- i__3, i__4, i__5;
-
- /* Builtin functions */
- integer f_rew(), s_wsfe(), e_wsfe(), do_fio();
-
- /* Local variables */
- static doublereal beta;
- static integer ldas, ldcs;
- static logical same;
- static doublereal bets;
- static logical tran, null;
- static char uplo[1];
- static integer i__, j, k, n;
- extern /* Subroutine */ int dmake_();
- static doublereal alpha;
- extern /* Subroutine */ int dmmch_();
- static logical isame[13];
- static integer nargs;
- static logical reset;
- static char trans[1];
- static logical upper;
- static char uplos[1];
- extern /* Subroutine */ void dprcn4_();
- static integer ia, ib, jc, ma, na, nc, ik, in, jj, lj, ks, ns;
- extern logical lderes_();
- static doublereal errmax;
- extern /* Subroutine */ int cdsyrk_();
- static char transs[1];
- static integer laa, lda, lcc, ldc;
- extern logical lde_();
- static doublereal als;
- static integer ict, icu;
- static doublereal err;
-
- /* Tests DSYRK. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Local Arrays .. */
- /* .. External Functions .. */
- /* .. External Subroutines .. */
- /* .. Intrinsic Functions .. */
- /* .. Scalars in Common .. */
- /* .. Common blocks .. */
- /* .. Data statements .. */
- /* Parameter adjustments */
- --idim;
- --alf;
- --bet;
- --g;
- --ct;
- --cs;
- --cc;
- c_dim1 = *nmax;
- c_offset = 1 + c_dim1 * 1;
- c__ -= c_offset;
- --bs;
- --bb;
- b_dim1 = *nmax;
- b_offset = 1 + b_dim1 * 1;
- b -= b_offset;
- --as;
- --aa;
- a_dim1 = *nmax;
- a_offset = 1 + a_dim1 * 1;
- a -= a_offset;
-
- /* Function Body */
- /* .. Executable Statements .. */
-
- nargs = 10;
- nc = 0;
- reset = TRUE_;
- errmax = 0.;
-
- i__1 = *nidim;
- for (in = 1; in <= i__1; ++in) {
- n = idim[in];
- /* Set LDC to 1 more than minimum value if room. */
- ldc = n;
- if (ldc < *nmax) {
- ++ldc;
- }
- /* Skip tests if not enough room. */
- if (ldc > *nmax) {
- goto L100;
- }
- lcc = ldc * n;
- null = n <= 0;
-
- i__2 = *nidim;
- for (ik = 1; ik <= i__2; ++ik) {
- k = idim[ik];
-
- for (ict = 1; ict <= 3; ++ict) {
- *(unsigned char *)trans = *(unsigned char *)&icht[ict - 1];
- tran = *(unsigned char *)trans == 'T' || *(unsigned char *)
- trans == 'C';
- if (tran) {
- ma = k;
- na = n;
- } else {
- ma = n;
- na = k;
- }
- /* Set LDA to 1 more than minimum value if room. */
- lda = ma;
- if (lda < *nmax) {
- ++lda;
- }
- /* Skip tests if not enough room. */
- if (lda > *nmax) {
- goto L80;
- }
- laa = lda * na;
-
- /* Generate the matrix A. */
-
- dmake_("GE", " ", " ", &ma, &na, &a[a_offset], nmax, &aa[1], &
- lda, &reset, &c_b104, (ftnlen)2, (ftnlen)1, (ftnlen)1)
- ;
-
- for (icu = 1; icu <= 2; ++icu) {
- *(unsigned char *)uplo = *(unsigned char *)&ichu[icu - 1];
- upper = *(unsigned char *)uplo == 'U';
-
- i__3 = *nalf;
- for (ia = 1; ia <= i__3; ++ia) {
- alpha = alf[ia];
-
- i__4 = *nbet;
- for (ib = 1; ib <= i__4; ++ib) {
- beta = bet[ib];
-
- /* Generate the matrix C. */
-
- dmake_("SY", uplo, " ", &n, &n, &c__[c_offset],
- nmax, &cc[1], &ldc, &reset, &c_b104, (
- ftnlen)2, (ftnlen)1, (ftnlen)1);
-
- ++nc;
-
- /* Save every datum before calling the subroutine. */
-
- *(unsigned char *)uplos = *(unsigned char *)uplo;
- *(unsigned char *)transs = *(unsigned char *)
- trans;
- ns = n;
- ks = k;
- als = alpha;
- i__5 = laa;
- for (i__ = 1; i__ <= i__5; ++i__) {
- as[i__] = aa[i__];
- /* L10: */
- }
- ldas = lda;
- bets = beta;
- i__5 = lcc;
- for (i__ = 1; i__ <= i__5; ++i__) {
- cs[i__] = cc[i__];
- /* L20: */
- }
- ldcs = ldc;
-
- /* Call the subroutine. */
-
- if (*trace) {
- dprcn4_(ntra, &nc, sname, iorder, uplo, trans,
- &n, &k, &alpha, &lda, &beta, &ldc, (
- ftnlen)12, (ftnlen)1, (ftnlen)1);
- }
- if (*rewi) {
- /* al__1.aerr = 0;
- al__1.aunit = *ntra;
- f_rew(&al__1);*/
- }
- cdsyrk_(iorder, uplo, trans, &n, &k, &alpha, &aa[
- 1], &lda, &beta, &cc[1], &ldc, (ftnlen)1,
- (ftnlen)1);
-
- /* Check if error-exit was taken incorrectly. */
-
- if (! infoc_1.ok) {
- printf("*** FATAL ERROR - ERROR-CALL MYEXIT TAKEN ON VALID CALL\n");
- *fatal = TRUE_;
- goto L120;
- }
-
- /* See what data changed inside subroutines. */
-
- isame[0] = *(unsigned char *)uplos == *(unsigned
- char *)uplo;
- isame[1] = *(unsigned char *)transs == *(unsigned
- char *)trans;
- isame[2] = ns == n;
- isame[3] = ks == k;
- isame[4] = als == alpha;
- isame[5] = lde_(&as[1], &aa[1], &laa);
- isame[6] = ldas == lda;
- isame[7] = bets == beta;
- if (null) {
- isame[8] = lde_(&cs[1], &cc[1], &lcc);
- } else {
- isame[8] = lderes_("SY", uplo, &n, &n, &cs[1],
- &cc[1], &ldc, (ftnlen)2, (ftnlen)1);
- }
- isame[9] = ldcs == ldc;
-
- /* If data was incorrectly changed, report and */
- /* return. */
-
- same = TRUE_;
- i__5 = nargs;
- for (i__ = 1; i__ <= i__5; ++i__) {
- same = same && isame[i__ - 1];
- if (! isame[i__ - 1]) {
- printf(" ******* FATAL ERROR - PARAMETER NUMBER %d WAS CHANGED INCORRECTLY *******\n",i__);
- }
- /* L30: */
- }
- if (! same) {
- *fatal = TRUE_;
- goto L120;
- }
-
- if (! null) {
-
- /* Check the result column by column. */
-
- jc = 1;
- i__5 = n;
- for (j = 1; j <= i__5; ++j) {
- if (upper) {
- jj = 1;
- lj = j;
- } else {
- jj = j;
- lj = n - j + 1;
- }
- if (tran) {
- dmmch_("T", "N", &lj, &c__1, &k, &
- alpha, &a[jj * a_dim1 + 1],
- nmax, &a[j * a_dim1 + 1],
- nmax, &beta, &c__[jj + j *
- c_dim1], nmax, &ct[1], &g[1],
- &cc[jc], &ldc, eps, &err,
- fatal, nout, &c_true, (ftnlen)
- 1, (ftnlen)1);
- } else {
- dmmch_("N", "T", &lj, &c__1, &k, &
- alpha, &a[jj + a_dim1], nmax,
- &a[j + a_dim1], nmax, &beta, &
- c__[jj + j * c_dim1], nmax, &
- ct[1], &g[1], &cc[jc], &ldc,
- eps, &err, fatal, nout, &
- c_true, (ftnlen)1, (ftnlen)1);
- }
- if (upper) {
- jc += ldc;
- } else {
- jc = jc + ldc + 1;
- }
- errmax = f2cmax(errmax,err);
- /* If got really bad answer, report and */
- /* return. */
- if (*fatal) {
- goto L110;
- }
- /* L40: */
- }
- }
-
- /* L50: */
- }
-
- /* L60: */
- }
-
- /* L70: */
- }
-
- L80:
- ;
- }
-
- /* L90: */
- }
-
- L100:
- ;
- }
-
- /* Report result. */
-
- if (errmax < *thresh) {
- if (*iorder == 0) {
- printf("%s PASSED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- if (*iorder == 1) {
- printf("%s PASSED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- } else {
- if (*iorder == 0) {
- printf("%s COMPLETED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- if (*iorder == 1) {
- printf("%s COMPLETED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- }
- goto L130;
-
- L110:
- if (n > 1) {
- printf(" THESE ARE THE RESULTS FOR COLUMN %d:\n",j);
- }
-
- L120:
- printf(" ******* %s FAILED ON CALL NUMBER:\n",sname);
- dprcn4_(nout, &nc, sname, iorder, uplo, trans, &n, &k, &alpha, &lda, &
- beta, &ldc, (ftnlen)12, (ftnlen)1, (ftnlen)1);
-
- L130:
- return 0;
-
- /* 9994 FORMAT( 1X, I6, ': ', A12,'(', 2( '''', A1, ''',' ), 2( I3, ',' ), */
- /* $ F4.1, ', A,', I3, ',', F4.1, ', C,', I3, ') .' ) */
-
- /* End of DCHK4. */
-
- } /* dchk4_ */
-
-
- /* Subroutine */ void dprcn4_(nout, nc, sname, iorder, uplo, transa, n, k,
- alpha, lda, beta, ldc, sname_len, uplo_len, transa_len)
- integer *nout, *nc;
- char *sname;
- integer *iorder;
- char *uplo, *transa;
- integer *n, *k;
- doublereal *alpha;
- integer *lda;
- doublereal *beta;
- integer *ldc;
- ftnlen sname_len;
- ftnlen uplo_len;
- ftnlen transa_len;
- {
- /* Builtin functions */
- integer s_wsfe(), do_fio(), e_wsfe();
-
- /* Local variables */
- static char ca[14], cu[14], crc[14];
-
- if (*(unsigned char *)uplo == 'U') {
- s_copy(cu, " CblasUpper", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cu, " CblasLower", (ftnlen)14, (ftnlen)14);
- }
- if (*(unsigned char *)transa == 'N') {
- s_copy(ca, " CblasNoTrans", (ftnlen)14, (ftnlen)14);
- } else if (*(unsigned char *)transa == 'T') {
- s_copy(ca, " CblasTrans", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(ca, "CblasConjTrans", (ftnlen)14, (ftnlen)14);
- }
- if (*iorder == 1) {
- s_copy(crc, " CblasRowMajor", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(crc, " CblasColMajor", (ftnlen)14, (ftnlen)14);
- }
- printf("%6d: %s %s %s %s\n",*nc,sname,crc,cu,ca);
- printf("( %d %d %4.1f A %d %4.1f C %d\n",*n,*k,*alpha,*lda,*beta,*ldc);
- } /* dprcn4_ */
-
-
- /* Subroutine */ int dchk5_(sname, eps, thresh, nout, ntra, trace, rewi,
- fatal, nidim, idim, nalf, alf, nbet, bet, nmax, ab, aa, as, bb, bs,
- c__, cc, cs, ct, g, w, iorder, sname_len)
- char *sname;
- doublereal *eps, *thresh;
- integer *nout, *ntra;
- logical *trace, *rewi, *fatal;
- integer *nidim, *idim, *nalf;
- doublereal *alf;
- integer *nbet;
- doublereal *bet;
- integer *nmax;
- doublereal *ab, *aa, *as, *bb, *bs, *c__, *cc, *cs, *ct, *g, *w;
- integer *iorder;
- ftnlen sname_len;
- {
- /* Initialized data */
-
- static char icht[3+1] = "NTC";
- static char ichu[2+1] = "UL";
-
- /* System generated locals */
- integer c_dim1, c_offset, i__1, i__2, i__3, i__4, i__5, i__6, i__7, i__8;
-
- /* Builtin functions */
- integer f_rew(), s_wsfe(), e_wsfe(), do_fio();
-
- /* Local variables */
- static integer jjab;
- static doublereal beta;
- static integer ldas, ldbs, ldcs;
- static logical same;
- static doublereal bets;
- static logical tran, null;
- static char uplo[1];
- static integer i__, j, k, n;
- extern /* Subroutine */ int dmake_();
- static doublereal alpha;
- extern /* Subroutine */ int dmmch_();
- static logical isame[13];
- static integer nargs;
- static logical reset;
- static char trans[1];
- static logical upper;
- static char uplos[1];
- extern /* Subroutine */ void dprcn5_();
- static integer ia, ib, jc, ma, na, nc, ik, in, jj, lj, ks, ns;
- extern logical lderes_();
- static doublereal errmax;
- static char transs[1];
- static integer laa, lbb, lda, lcc, ldb, ldc;
- extern logical lde_();
- extern /* Subroutine */ int cdsyr2k_();
- static doublereal als;
- static integer ict, icu;
- static doublereal err;
-
- /* Tests DSYR2K. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Local Arrays .. */
- /* .. External Functions .. */
- /* .. External Subroutines .. */
- /* .. Intrinsic Functions .. */
- /* .. Scalars in Common .. */
- /* .. Common blocks .. */
- /* .. Data statements .. */
- /* Parameter adjustments */
- --idim;
- --alf;
- --bet;
- --w;
- --g;
- --ct;
- --cs;
- --cc;
- c_dim1 = *nmax;
- c_offset = 1 + c_dim1 * 1;
- c__ -= c_offset;
- --bs;
- --bb;
- --as;
- --aa;
- --ab;
-
- /* Function Body */
- /* .. Executable Statements .. */
-
- nargs = 12;
- nc = 0;
- reset = TRUE_;
- errmax = 0.;
-
- i__1 = *nidim;
- for (in = 1; in <= i__1; ++in) {
- n = idim[in];
- /* Set LDC to 1 more than minimum value if room. */
- ldc = n;
- if (ldc < *nmax) {
- ++ldc;
- }
- /* Skip tests if not enough room. */
- if (ldc > *nmax) {
- goto L130;
- }
- lcc = ldc * n;
- null = n <= 0;
-
- i__2 = *nidim;
- for (ik = 1; ik <= i__2; ++ik) {
- k = idim[ik];
-
- for (ict = 1; ict <= 3; ++ict) {
- *(unsigned char *)trans = *(unsigned char *)&icht[ict - 1];
- tran = *(unsigned char *)trans == 'T' || *(unsigned char *)
- trans == 'C';
- if (tran) {
- ma = k;
- na = n;
- } else {
- ma = n;
- na = k;
- }
- /* Set LDA to 1 more than minimum value if room. */
- lda = ma;
- if (lda < *nmax) {
- ++lda;
- }
- /* Skip tests if not enough room. */
- if (lda > *nmax) {
- goto L110;
- }
- laa = lda * na;
-
- /* Generate the matrix A. */
-
- if (tran) {
- i__3 = *nmax << 1;
- dmake_("GE", " ", " ", &ma, &na, &ab[1], &i__3, &aa[1], &
- lda, &reset, &c_b104, (ftnlen)2, (ftnlen)1, (
- ftnlen)1);
- } else {
- dmake_("GE", " ", " ", &ma, &na, &ab[1], nmax, &aa[1], &
- lda, &reset, &c_b104, (ftnlen)2, (ftnlen)1, (
- ftnlen)1);
- }
-
- /* Generate the matrix B. */
-
- ldb = lda;
- lbb = laa;
- if (tran) {
- i__3 = *nmax << 1;
- dmake_("GE", " ", " ", &ma, &na, &ab[k + 1], &i__3, &bb[1]
- , &ldb, &reset, &c_b104, (ftnlen)2, (ftnlen)1, (
- ftnlen)1);
- } else {
- dmake_("GE", " ", " ", &ma, &na, &ab[k * *nmax + 1], nmax,
- &bb[1], &ldb, &reset, &c_b104, (ftnlen)2, (
- ftnlen)1, (ftnlen)1);
- }
-
- for (icu = 1; icu <= 2; ++icu) {
- *(unsigned char *)uplo = *(unsigned char *)&ichu[icu - 1];
- upper = *(unsigned char *)uplo == 'U';
-
- i__3 = *nalf;
- for (ia = 1; ia <= i__3; ++ia) {
- alpha = alf[ia];
-
- i__4 = *nbet;
- for (ib = 1; ib <= i__4; ++ib) {
- beta = bet[ib];
-
- /* Generate the matrix C. */
-
- dmake_("SY", uplo, " ", &n, &n, &c__[c_offset],
- nmax, &cc[1], &ldc, &reset, &c_b104, (
- ftnlen)2, (ftnlen)1, (ftnlen)1);
-
- ++nc;
-
- /* Save every datum before calling the subroutine. */
-
- *(unsigned char *)uplos = *(unsigned char *)uplo;
- *(unsigned char *)transs = *(unsigned char *)
- trans;
- ns = n;
- ks = k;
- als = alpha;
- i__5 = laa;
- for (i__ = 1; i__ <= i__5; ++i__) {
- as[i__] = aa[i__];
- /* L10: */
- }
- ldas = lda;
- i__5 = lbb;
- for (i__ = 1; i__ <= i__5; ++i__) {
- bs[i__] = bb[i__];
- /* L20: */
- }
- ldbs = ldb;
- bets = beta;
- i__5 = lcc;
- for (i__ = 1; i__ <= i__5; ++i__) {
- cs[i__] = cc[i__];
- /* L30: */
- }
- ldcs = ldc;
-
- /* Call the subroutine. */
-
- if (*trace) {
- dprcn5_(ntra, &nc, sname, iorder, uplo, trans,
- &n, &k, &alpha, &lda, &ldb, &beta, &
- ldc, (ftnlen)12, (ftnlen)1, (ftnlen)1)
- ;
- }
- if (*rewi) {
- /* al__1.aerr = 0;
- al__1.aunit = *ntra;
- f_rew(&al__1);*/
- }
- cdsyr2k_(iorder, uplo, trans, &n, &k, &alpha, &aa[
- 1], &lda, &bb[1], &ldb, &beta, &cc[1], &
- ldc, (ftnlen)1, (ftnlen)1);
-
- /* Check if error-exit was taken incorrectly. */
-
- if (! infoc_1.ok) {
- printf("*** FATAL ERROR - ERROR-CALL MYEXIT TAKEN ON VALID CALL\n");
- *fatal = TRUE_;
- goto L150;
- }
-
- /* See what data changed inside subroutines. */
-
- isame[0] = *(unsigned char *)uplos == *(unsigned
- char *)uplo;
- isame[1] = *(unsigned char *)transs == *(unsigned
- char *)trans;
- isame[2] = ns == n;
- isame[3] = ks == k;
- isame[4] = als == alpha;
- isame[5] = lde_(&as[1], &aa[1], &laa);
- isame[6] = ldas == lda;
- isame[7] = lde_(&bs[1], &bb[1], &lbb);
- isame[8] = ldbs == ldb;
- isame[9] = bets == beta;
- if (null) {
- isame[10] = lde_(&cs[1], &cc[1], &lcc);
- } else {
- isame[10] = lderes_("SY", uplo, &n, &n, &cs[1]
- , &cc[1], &ldc, (ftnlen)2, (ftnlen)1);
- }
- isame[11] = ldcs == ldc;
-
- /* If data was incorrectly changed, report and */
- /* return. */
-
- same = TRUE_;
- i__5 = nargs;
- for (i__ = 1; i__ <= i__5; ++i__) {
- same = same && isame[i__ - 1];
- if (! isame[i__ - 1]) {
- printf(" ******* FATAL ERROR - PARAMETER NUMBER %d WAS CHANGED INCORRECTLY *******\n",i__);
- }
- /* L40: */
- }
- if (! same) {
- *fatal = TRUE_;
- goto L150;
- }
-
- if (! null) {
-
- /* Check the result column by column. */
-
- jjab = 1;
- jc = 1;
- i__5 = n;
- for (j = 1; j <= i__5; ++j) {
- if (upper) {
- jj = 1;
- lj = j;
- } else {
- jj = j;
- lj = n - j + 1;
- }
- if (tran) {
- i__6 = k;
- for (i__ = 1; i__ <= i__6; ++i__) {
- w[i__] = ab[((j - 1) << 1) * *nmax
- + k + i__];
- w[k + i__] = ab[((j - 1) << 1) * *
- nmax + i__];
- /* L50: */
- }
- i__6 = k << 1;
- i__7 = *nmax << 1;
- i__8 = *nmax << 1;
- dmmch_("T", "N", &lj, &c__1, &i__6, &
- alpha, &ab[jjab], &i__7, &w[1]
- , &i__8, &beta, &c__[jj + j *
- c_dim1], nmax, &ct[1], &g[1],
- &cc[jc], &ldc, eps, &err,
- fatal, nout, &c_true, (ftnlen)
- 1, (ftnlen)1);
- } else {
- i__6 = k;
- for (i__ = 1; i__ <= i__6; ++i__) {
- w[i__] = ab[(k + i__ - 1) * *nmax
- + j];
- w[k + i__] = ab[(i__ - 1) * *nmax
- + j];
- /* L60: */
- }
- i__6 = k << 1;
- i__7 = *nmax << 1;
- dmmch_("N", "N", &lj, &c__1, &i__6, &
- alpha, &ab[jj], nmax, &w[1], &
- i__7, &beta, &c__[jj + j *
- c_dim1], nmax, &ct[1], &g[1],
- &cc[jc], &ldc, eps, &err,
- fatal, nout, &c_true, (ftnlen)
- 1, (ftnlen)1);
- }
- if (upper) {
- jc += ldc;
- } else {
- jc = jc + ldc + 1;
- if (tran) {
- jjab += *nmax << 1;
- }
- }
- errmax = f2cmax(errmax,err);
- /* If got really bad answer, report and */
- /* return. */
- if (*fatal) {
- goto L140;
- }
- /* L70: */
- }
- }
-
- /* L80: */
- }
-
- /* L90: */
- }
-
- /* L100: */
- }
-
- L110:
- ;
- }
-
- /* L120: */
- }
-
- L130:
- ;
- }
-
- /* Report result. */
-
- if (errmax < *thresh) {
- if (*iorder == 0) {
- printf("%s PASSED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- if (*iorder == 1) {
- printf("%s PASSED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)\n",sname,nc);
- }
- } else {
- if (*iorder == 0) {
- printf("%s COMPLETED THE COLUMN-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- if (*iorder == 1) {
- printf("%s COMPLETED THE ROW-MAJOR COMPUTATIONAL TESTS (%d CALLS)/n",sname,nc);
- printf("***** BUT WITH MAXIMUM TEST RATIO %8.2f - SUSPECT *******/n",errmax);
- }
- }
- goto L160;
-
- L140:
- if (n > 1) {
- printf(" THESE ARE THE RESULTS FOR COLUMN %d:\n",j);
- }
-
- L150:
- printf(" ******* %s FAILED ON CALL NUMBER:\n",sname);
- dprcn5_(nout, &nc, sname, iorder, uplo, trans, &n, &k, &alpha, &lda, &ldb,
- &beta, &ldc, (ftnlen)12, (ftnlen)1, (ftnlen)1);
-
- L160:
- return 0;
-
- /* 9994 FORMAT( 1X, I6, ': ', A12,'(', 2( '''', A1, ''',' ), 2( I3, ',' ), */
- /* $ F4.1, ', A,', I3, ', B,', I3, ',', F4.1, ', C,', I3, ') ', */
- /* $ ' .' ) */
-
- /* End of DCHK5. */
-
- } /* dchk5_ */
-
-
- /* Subroutine */ void dprcn5_(nout, nc, sname, iorder, uplo, transa, n, k,
- alpha, lda, ldb, beta, ldc, sname_len, uplo_len, transa_len)
- integer *nout, *nc;
- char *sname;
- integer *iorder;
- char *uplo, *transa;
- integer *n, *k;
- doublereal *alpha;
- integer *lda, *ldb;
- doublereal *beta;
- integer *ldc;
- ftnlen sname_len;
- ftnlen uplo_len;
- ftnlen transa_len;
- {
- /* Builtin functions */
- integer s_wsfe(), do_fio(), e_wsfe();
-
- /* Local variables */
- static char ca[14], cu[14], crc[14];
-
- if (*(unsigned char *)uplo == 'U') {
- s_copy(cu, " CblasUpper", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(cu, " CblasLower", (ftnlen)14, (ftnlen)14);
- }
- if (*(unsigned char *)transa == 'N') {
- s_copy(ca, " CblasNoTrans", (ftnlen)14, (ftnlen)14);
- } else if (*(unsigned char *)transa == 'T') {
- s_copy(ca, " CblasTrans", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(ca, "CblasConjTrans", (ftnlen)14, (ftnlen)14);
- }
- if (*iorder == 1) {
- s_copy(crc, " CblasRowMajor", (ftnlen)14, (ftnlen)14);
- } else {
- s_copy(crc, " CblasColMajor", (ftnlen)14, (ftnlen)14);
- }
- printf("%6d: %s %s %s %s\n",*nc,sname,crc,cu,ca);
- printf("%d %d %4.1f , A, %d, B, %d, %4.1f , C, %d.\n",*n,*k,*alpha,*lda,*ldb,*beta,*ldc);
- } /* dprcn5_ */
-
-
- /* Subroutine */ int dmake_(type__, uplo, diag, m, n, a, nmax, aa, lda, reset,
- transl, type_len, uplo_len, diag_len)
- char *type__, *uplo, *diag;
- integer *m, *n;
- doublereal *a;
- integer *nmax;
- doublereal *aa;
- integer *lda;
- logical *reset;
- doublereal *transl;
- ftnlen type_len;
- ftnlen uplo_len;
- ftnlen diag_len;
- {
- /* System generated locals */
- integer a_dim1, a_offset, i__1, i__2;
-
- /* Local variables */
- extern doublereal dbeg_();
- static integer ibeg, iend;
- static logical unit;
- static integer i__, j;
- static logical lower, upper, gen, tri, sym;
-
-
- /* Generates values for an M by N matrix A. */
- /* Stores the values in the array AA in the data structure required */
- /* by the routine, with unwanted elements set to rogue value. */
-
- /* TYPE is 'GE', 'SY' or 'TR'. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. External Functions .. */
- /* .. Executable Statements .. */
- /* Parameter adjustments */
- a_dim1 = *nmax;
- a_offset = 1 + a_dim1 * 1;
- a -= a_offset;
- --aa;
-
- /* Function Body */
- gen = s_cmp(type__, "GE", (ftnlen)2, (ftnlen)2) == 0;
- sym = s_cmp(type__, "SY", (ftnlen)2, (ftnlen)2) == 0;
- tri = s_cmp(type__, "TR", (ftnlen)2, (ftnlen)2) == 0;
- upper = (sym || tri) && *(unsigned char *)uplo == 'U';
- lower = (sym || tri) && *(unsigned char *)uplo == 'L';
- unit = tri && *(unsigned char *)diag == 'U';
-
- /* Generate data in array A. */
-
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
- i__2 = *m;
- for (i__ = 1; i__ <= i__2; ++i__) {
- if (gen || (upper && i__ <= j) || (lower && i__ >= j)) {
- a[i__ + j * a_dim1] = dbeg_(reset) + *transl;
- if (i__ != j) {
- /* Set some elements to zero */
- if (*n > 3 && j == *n / 2) {
- a[i__ + j * a_dim1] = 0.;
- }
- if (sym) {
- a[j + i__ * a_dim1] = a[i__ + j * a_dim1];
- } else if (tri) {
- a[j + i__ * a_dim1] = 0.;
- }
- }
- }
- /* L10: */
- }
- if (tri) {
- a[j + j * a_dim1] += 1.;
- }
- if (unit) {
- a[j + j * a_dim1] = 1.;
- }
- /* L20: */
- }
-
- /* Store elements in array AS in data structure required by routine. */
-
- if (s_cmp(type__, "GE", (ftnlen)2, (ftnlen)2) == 0) {
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
- i__2 = *m;
- for (i__ = 1; i__ <= i__2; ++i__) {
- aa[i__ + (j - 1) * *lda] = a[i__ + j * a_dim1];
- /* L30: */
- }
- i__2 = *lda;
- for (i__ = *m + 1; i__ <= i__2; ++i__) {
- aa[i__ + (j - 1) * *lda] = -1e10;
- /* L40: */
- }
- /* L50: */
- }
- } else if (s_cmp(type__, "SY", (ftnlen)2, (ftnlen)2) == 0 || s_cmp(type__,
- "TR", (ftnlen)2, (ftnlen)2) == 0) {
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
- if (upper) {
- ibeg = 1;
- if (unit) {
- iend = j - 1;
- } else {
- iend = j;
- }
- } else {
- if (unit) {
- ibeg = j + 1;
- } else {
- ibeg = j;
- }
- iend = *n;
- }
- i__2 = ibeg - 1;
- for (i__ = 1; i__ <= i__2; ++i__) {
- aa[i__ + (j - 1) * *lda] = -1e10;
- /* L60: */
- }
- i__2 = iend;
- for (i__ = ibeg; i__ <= i__2; ++i__) {
- aa[i__ + (j - 1) * *lda] = a[i__ + j * a_dim1];
- /* L70: */
- }
- i__2 = *lda;
- for (i__ = iend + 1; i__ <= i__2; ++i__) {
- aa[i__ + (j - 1) * *lda] = -1e10;
- /* L80: */
- }
- /* L90: */
- }
- }
- return 0;
-
- /* End of DMAKE. */
-
- } /* dmake_ */
-
- /* Subroutine */ int dmmch_(transa, transb, m, n, kk, alpha, a, lda, b, ldb,
- beta, c__, ldc, ct, g, cc, ldcc, eps, err, fatal, nout, mv,
- transa_len, transb_len)
- char *transa, *transb;
- integer *m, *n, *kk;
- doublereal *alpha, *a;
- integer *lda;
- doublereal *b;
- integer *ldb;
- doublereal *beta, *c__;
- integer *ldc;
- doublereal *ct, *g, *cc;
- integer *ldcc;
- doublereal *eps, *err;
- logical *fatal;
- integer *nout;
- logical *mv;
- ftnlen transa_len;
- ftnlen transb_len;
- {
- /* System generated locals */
- integer a_dim1, a_offset, b_dim1, b_offset, c_dim1, c_offset, cc_dim1,
- cc_offset, i__1, i__2, i__3;
- doublereal d__1, d__2;
-
- /* Builtin functions */
- double sqrt();
- integer s_wsfe(), e_wsfe(), do_fio();
-
- /* Local variables */
- static doublereal erri;
- static integer i__, j, k;
- static logical trana, tranb;
-
- /* Checks the results of the computational tests. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Parameters .. */
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Intrinsic Functions .. */
- /* .. Executable Statements .. */
- /* Parameter adjustments */
- a_dim1 = *lda;
- a_offset = 1 + a_dim1 * 1;
- a -= a_offset;
- b_dim1 = *ldb;
- b_offset = 1 + b_dim1 * 1;
- b -= b_offset;
- c_dim1 = *ldc;
- c_offset = 1 + c_dim1 * 1;
- c__ -= c_offset;
- --ct;
- --g;
- cc_dim1 = *ldcc;
- cc_offset = 1 + cc_dim1 * 1;
- cc -= cc_offset;
-
- /* Function Body */
- trana = *(unsigned char *)transa == 'T' || *(unsigned char *)transa ==
- 'C';
- tranb = *(unsigned char *)transb == 'T' || *(unsigned char *)transb ==
- 'C';
-
- /* Compute expected result, one column at a time, in CT using data */
- /* in A, B and C. */
- /* Compute gauges in G. */
-
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
-
- i__2 = *m;
- for (i__ = 1; i__ <= i__2; ++i__) {
- ct[i__] = 0.;
- g[i__] = 0.;
- /* L10: */
- }
- if (! trana && ! tranb) {
- i__2 = *kk;
- for (k = 1; k <= i__2; ++k) {
- i__3 = *m;
- for (i__ = 1; i__ <= i__3; ++i__) {
- ct[i__] += a[i__ + k * a_dim1] * b[k + j * b_dim1];
- g[i__] += (d__1 = a[i__ + k * a_dim1], abs(d__1)) * (d__2
- = b[k + j * b_dim1], abs(d__2));
- /* L20: */
- }
- /* L30: */
- }
- } else if (trana && ! tranb) {
- i__2 = *kk;
- for (k = 1; k <= i__2; ++k) {
- i__3 = *m;
- for (i__ = 1; i__ <= i__3; ++i__) {
- ct[i__] += a[k + i__ * a_dim1] * b[k + j * b_dim1];
- g[i__] += (d__1 = a[k + i__ * a_dim1], abs(d__1)) * (d__2
- = b[k + j * b_dim1], abs(d__2));
- /* L40: */
- }
- /* L50: */
- }
- } else if (! trana && tranb) {
- i__2 = *kk;
- for (k = 1; k <= i__2; ++k) {
- i__3 = *m;
- for (i__ = 1; i__ <= i__3; ++i__) {
- ct[i__] += a[i__ + k * a_dim1] * b[j + k * b_dim1];
- g[i__] += (d__1 = a[i__ + k * a_dim1], abs(d__1)) * (d__2
- = b[j + k * b_dim1], abs(d__2));
- /* L60: */
- }
- /* L70: */
- }
- } else if (trana && tranb) {
- i__2 = *kk;
- for (k = 1; k <= i__2; ++k) {
- i__3 = *m;
- for (i__ = 1; i__ <= i__3; ++i__) {
- ct[i__] += a[k + i__ * a_dim1] * b[j + k * b_dim1];
- g[i__] += (d__1 = a[k + i__ * a_dim1], abs(d__1)) * (d__2
- = b[j + k * b_dim1], abs(d__2));
- /* L80: */
- }
- /* L90: */
- }
- }
- i__2 = *m;
- for (i__ = 1; i__ <= i__2; ++i__) {
- ct[i__] = *alpha * ct[i__] + *beta * c__[i__ + j * c_dim1];
- g[i__] = abs(*alpha) * g[i__] + abs(*beta) * (d__1 = c__[i__ + j *
- c_dim1], abs(d__1));
- /* L100: */
- }
-
- /* Compute the error ratio for this result. */
-
- *err = 0.;
- i__2 = *m;
- for (i__ = 1; i__ <= i__2; ++i__) {
- erri = (d__1 = ct[i__] - cc[i__ + j * cc_dim1], abs(d__1)) / *eps;
- if (g[i__] != 0.) {
- erri /= g[i__];
- }
- *err = f2cmax(*err,erri);
- if (*err * sqrt(*eps) >= 1.) {
- goto L130;
- }
- /* L110: */
- }
-
- /* L120: */
- }
-
- /* If the loop completes, all results are at least half accurate. */
- goto L150;
-
- /* Report fatal error. */
-
- L130:
- *fatal = TRUE_;
- printf(" ******* FATAL ERROR - COMPUTED RESULT IS LESS THAN HALF ACCURATE *******\n");
- printf(" EXPECTED RESULT COMPUTED RESULT\n");
- i__1 = *m;
- for (i__ = 1; i__ <= i__1; ++i__) {
- if (*mv) {
- printf("%7d %15.6g %15.6g\n",i__,ct[i__],cc[i__+j*cc_dim1]);
- } else {
- printf("%7d %15.6g %15.6g\n",i__,cc[i__+j*cc_dim1],ct[i__]);
- }
- /* L140: */
- }
- if (*n > 1) {
- printf(" THESE ARE THE RESULTS FOR COLUMN %d\n",j);
- }
-
- L150:
- return 0;
-
-
- /* End of DMMCH. */
-
- } /* dmmch_ */
-
- logical lde_(ri, rj, lr)
- doublereal *ri, *rj;
- integer *lr;
- {
- /* System generated locals */
- integer i__1;
- logical ret_val;
-
- /* Local variables */
- static integer i__;
-
-
- /* Tests if two arrays are identical. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Executable Statements .. */
- /* Parameter adjustments */
- --rj;
- --ri;
-
- /* Function Body */
- i__1 = *lr;
- for (i__ = 1; i__ <= i__1; ++i__) {
- if (ri[i__] != rj[i__]) {
- goto L20;
- }
- /* L10: */
- }
- ret_val = TRUE_;
- goto L30;
- L20:
- ret_val = FALSE_;
- L30:
- return ret_val;
-
- /* End of LDE. */
-
- } /* lde_ */
-
- logical lderes_(type__, uplo, m, n, aa, as, lda, type_len, uplo_len)
- char *type__, *uplo;
- integer *m, *n;
- doublereal *aa, *as;
- integer *lda;
- ftnlen type_len;
- ftnlen uplo_len;
- {
- /* System generated locals */
- integer aa_dim1, aa_offset, as_dim1, as_offset, i__1, i__2;
- logical ret_val;
-
- /* Local variables */
- static integer ibeg, iend, i__, j;
- static logical upper;
-
-
- /* Tests if selected elements in two arrays are equal. */
-
- /* TYPE is 'GE' or 'SY'. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Scalar Arguments .. */
- /* .. Array Arguments .. */
- /* .. Local Scalars .. */
- /* .. Executable Statements .. */
- /* Parameter adjustments */
- as_dim1 = *lda;
- as_offset = 1 + as_dim1 * 1;
- as -= as_offset;
- aa_dim1 = *lda;
- aa_offset = 1 + aa_dim1 * 1;
- aa -= aa_offset;
-
- /* Function Body */
- upper = *(unsigned char *)uplo == 'U';
- if (s_cmp(type__, "GE", (ftnlen)2, (ftnlen)2) == 0) {
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
- i__2 = *lda;
- for (i__ = *m + 1; i__ <= i__2; ++i__) {
- if (aa[i__ + j * aa_dim1] != as[i__ + j * as_dim1]) {
- goto L70;
- }
- /* L10: */
- }
- /* L20: */
- }
- } else if (s_cmp(type__, "SY", (ftnlen)2, (ftnlen)2) == 0) {
- i__1 = *n;
- for (j = 1; j <= i__1; ++j) {
- if (upper) {
- ibeg = 1;
- iend = j;
- } else {
- ibeg = j;
- iend = *n;
- }
- i__2 = ibeg - 1;
- for (i__ = 1; i__ <= i__2; ++i__) {
- if (aa[i__ + j * aa_dim1] != as[i__ + j * as_dim1]) {
- goto L70;
- }
- /* L30: */
- }
- i__2 = *lda;
- for (i__ = iend + 1; i__ <= i__2; ++i__) {
- if (aa[i__ + j * aa_dim1] != as[i__ + j * as_dim1]) {
- goto L70;
- }
- /* L40: */
- }
- /* L50: */
- }
- }
-
- /* 60 CONTINUE */
- ret_val = TRUE_;
- goto L80;
- L70:
- ret_val = FALSE_;
- L80:
- return ret_val;
-
- /* End of LDERES. */
-
- } /* lderes_ */
-
- doublereal dbeg_(reset)
- logical *reset;
- {
- /* System generated locals */
- doublereal ret_val;
-
- /* Local variables */
- static integer i__, ic, mi;
-
-
- /* Generates random numbers uniformly distributed between -0.5 and 0.5. */
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Scalar Arguments .. */
- /* .. Local Scalars .. */
- /* .. Save statement .. */
- /* .. Executable Statements .. */
- if (*reset) {
- /* Initialize local variables. */
- mi = 891;
- i__ = 7;
- ic = 0;
- *reset = FALSE_;
- }
-
- /* The sequence of values of I is bounded between 1 and 999. */
- /* If initial I = 1,2,3,6,7 or 9, the period will be 50. */
- /* If initial I = 4 or 8, the period will be 25. */
- /* If initial I = 5, the period will be 10. */
- /* IC is used to break up the period by skipping 1 value of I in 6. */
-
- ++ic;
- L10:
- i__ *= mi;
- i__ -= i__ / 1000 * 1000;
- if (ic >= 5) {
- ic = 0;
- goto L10;
- }
- ret_val = (i__ - 500) / 1001.;
- return ret_val;
-
- /* End of DBEG. */
-
- } /* dbeg_ */
-
- doublereal ddiff_(x, y)
- doublereal *x, *y;
- {
- /* System generated locals */
- doublereal ret_val;
-
-
- /* Auxiliary routine for test program for Level 3 Blas. */
-
- /* -- Written on 8-February-1989. */
- /* Jack Dongarra, Argonne National Laboratory. */
- /* Iain Duff, AERE Harwell. */
- /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
- /* Sven Hammarling, Numerical Algorithms Group Ltd. */
-
- /* .. Scalar Arguments .. */
- /* .. Executable Statements .. */
- ret_val = *x - *y;
- return ret_val;
-
- /* End of DDIFF. */
-
- } /* ddiff_ */
-
- /* Main program alias */ /*int dblat3_ () { MAIN__ (); }*/
|