| 123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167116811691170117111721173117411751176117711781179118011811182118311841185118611871188118911901191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260126112621263126412651266126712681269127012711272127312741275127612771278127912801281128212831284128512861287128812891290129112921293129412951296129712981299130013011302130313041305130613071308130913101311131213131314131513161317131813191320132113221323132413251326132713281329133013311332133313341335133613371338133913401341134213431344134513461347134813491350135113521353135413551356135713581359136013611362136313641365136613671368136913701371137213731374137513761377137813791380138113821383138413851386138713881389139013911392139313941395139613971398139914001401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470147114721473147414751476147714781479148014811482148314841485148614871488148914901491149214931494149514961497149814991500150115021503150415051506150715081509151015111512151315141515151615171518151915201521152215231524152515261527152815291530153115321533153415351536153715381539154015411542154315441545154615471548154915501551155215531554155515561557155815591560156115621563156415651566156715681569157015711572157315741575157615771578157915801581158215831584158515861587158815891590159115921593159415951596159715981599160016011602160316041605160616071608160916101611161216131614161516161617161816191620162116221623162416251626162716281629163016311632163316341635163616371638163916401641164216431644164516461647164816491650165116521653165416551656165716581659166016611662166316641665166616671668166916701671167216731674167516761677167816791680168116821683168416851686168716881689169016911692169316941695169616971698169917001701170217031704170517061707170817091710171117121713171417151716171717181719172017211722172317241725172617271728172917301731173217331734173517361737173817391740174117421743174417451746174717481749175017511752175317541755175617571758175917601761176217631764176517661767176817691770177117721773177417751776177717781779178017811782178317841785178617871788178917901791179217931794179517961797179817991800180118021803180418051806180718081809181018111812181318141815181618171818181918201821182218231824182518261827182818291830183118321833183418351836183718381839184018411842184318441845184618471848184918501851185218531854185518561857185818591860186118621863186418651866186718681869187018711872187318741875187618771878187918801881188218831884188518861887188818891890189118921893189418951896189718981899190019011902190319041905190619071908190919101911191219131914191519161917191819191920192119221923192419251926192719281929193019311932193319341935193619371938193919401941194219431944194519461947194819491950195119521953195419551956195719581959196019611962196319641965196619671968196919701971197219731974197519761977197819791980198119821983198419851986198719881989199019911992199319941995199619971998199920002001200220032004200520062007200820092010201120122013201420152016201720182019202020212022202320242025202620272028202920302031203220332034203520362037203820392040204120422043204420452046204720482049205020512052205320542055205620572058205920602061206220632064206520662067206820692070207120722073207420752076207720782079208020812082208320842085208620872088208920902091209220932094209520962097209820992100210121022103210421052106210721082109211021112112211321142115211621172118211921202121212221232124212521262127212821292130213121322133213421352136213721382139214021412142214321442145214621472148214921502151215221532154215521562157215821592160216121622163216421652166216721682169217021712172217321742175217621772178217921802181218221832184218521862187218821892190219121922193219421952196219721982199220022012202220322042205220622072208220922102211221222132214221522162217221822192220222122222223222422252226222722282229223022312232223322342235223622372238223922402241224222432244224522462247224822492250225122522253225422552256225722582259226022612262226322642265226622672268226922702271227222732274227522762277227822792280228122822283228422852286228722882289229022912292229322942295229622972298229923002301230223032304230523062307230823092310231123122313231423152316231723182319232023212322232323242325232623272328232923302331233223332334233523362337233823392340234123422343234423452346234723482349235023512352235323542355235623572358235923602361236223632364236523662367236823692370237123722373237423752376237723782379238023812382238323842385238623872388238923902391239223932394239523962397239823992400240124022403240424052406240724082409241024112412241324142415241624172418241924202421242224232424242524262427242824292430243124322433243424352436243724382439244024412442244324442445244624472448244924502451245224532454245524562457245824592460246124622463246424652466246724682469247024712472247324742475247624772478247924802481248224832484248524862487248824892490249124922493249424952496249724982499250025012502250325042505250625072508250925102511251225132514251525162517251825192520252125222523252425252526252725282529253025312532253325342535253625372538253925402541254225432544254525462547254825492550255125522553255425552556255725582559256025612562256325642565256625672568256925702571257225732574257525762577257825792580258125822583258425852586258725882589259025912592259325942595259625972598259926002601260226032604260526062607260826092610261126122613261426152616261726182619262026212622262326242625262626272628262926302631263226332634263526362637263826392640264126422643264426452646264726482649265026512652265326542655265626572658265926602661266226632664266526662667266826692670267126722673267426752676267726782679268026812682268326842685268626872688268926902691269226932694269526962697269826992700270127022703270427052706270727082709271027112712271327142715271627172718271927202721272227232724272527262727272827292730273127322733273427352736273727382739274027412742274327442745274627472748274927502751275227532754275527562757275827592760276127622763276427652766276727682769277027712772277327742775277627772778277927802781278227832784278527862787278827892790279127922793279427952796279727982799280028012802280328042805280628072808280928102811281228132814281528162817281828192820282128222823282428252826282728282829283028312832283328342835283628372838283928402841284228432844284528462847284828492850285128522853285428552856285728582859286028612862286328642865286628672868286928702871287228732874287528762877287828792880288128822883288428852886288728882889289028912892289328942895289628972898289929002901290229032904290529062907290829092910291129122913291429152916291729182919292029212922292329242925292629272928292929302931293229332934293529362937293829392940294129422943294429452946294729482949295029512952295329542955295629572958295929602961296229632964296529662967296829692970297129722973297429752976297729782979298029812982298329842985298629872988298929902991299229932994299529962997299829993000300130023003300430053006300730083009301030113012301330143015301630173018301930203021302230233024302530263027302830293030303130323033303430353036303730383039304030413042304330443045304630473048304930503051305230533054305530563057305830593060306130623063306430653066306730683069307030713072307330743075307630773078307930803081308230833084308530863087308830893090309130923093309430953096309730983099310031013102310331043105310631073108310931103111311231133114311531163117311831193120312131223123312431253126312731283129313031313132313331343135313631373138313931403141314231433144314531463147314831493150315131523153315431553156315731583159316031613162316331643165316631673168316931703171317231733174317531763177317831793180318131823183318431853186318731883189319031913192319331943195319631973198319932003201320232033204320532063207320832093210321132123213321432153216321732183219322032213222322332243225322632273228322932303231323232333234323532363237323832393240324132423243324432453246324732483249325032513252325332543255325632573258325932603261 |
- /*******************************************************************************
- * Copyright 2016-2022 Intel Corporation.
- *
- * This software and the related documents are Intel copyrighted materials, and
- * your use of them is governed by the express license under which they were
- * provided to you (License). Unless the License provides otherwise, you may not
- * use, modify, copy, publish, distribute, disclose or transmit this software or
- * the related documents without Intel's prior written permission.
- *
- * This software and the related documents are provided as is, with no express
- * or implied warranties, other than those that are expressly stated in the
- * License.
- *******************************************************************************/
- /*
- ! Content:
- ! Intel(R) oneAPI Math Kernel Library (oneMKL) C functions that can be inlined
- !******************************************************************************/
- #ifndef MKL_DC_UNSAFE
- #define MKL_DC_UNSAFE 0
- #endif
- #undef mkl_dc_getrf
- #undef mkl_dc_lapacke_getrf_convert
- #undef mkl_dc_getrf_1x1
- #undef mkl_dc_getrf_2x2
- #undef mkl_dc_getrf_3x3
- #undef mkl_dc_getrf_4x4
- #undef mkl_dc_getrf_5x5
- #undef mkl_dc_getrf_mxn
- #undef mkl_dc_getrf_generic
- #undef mkl_dc_getrfnp
- #undef mkl_dc_lapacke_getrfnp_convert
- #undef mkl_dc_getrfnp_1x1
- #undef mkl_dc_getrfnp_2x2
- #undef mkl_dc_getrfnp_3x3
- #undef mkl_dc_getrfnp_4x4
- #undef mkl_dc_getrfnp_5x5
- #undef mkl_dc_getrfnp_mxn
- #undef mkl_dc_getrfnp_generic
- #undef mkl_dc_getri
- #undef mkl_dc_lapacke_getri_convert
- #undef mkl_dc_getri_1x1
- #undef mkl_dc_getri_2x2
- #undef mkl_dc_getri_3x3
- #undef mkl_dc_getri_4x4
- #undef mkl_dc_getri_5x5
- #undef mkl_dc_getri_6x6
- #undef mkl_dc_getri_n
- #undef mkl_dc_getri_inverse_upper
- #undef mkl_dc_getri_solve_lower
- #undef mkl_dc_getri_generic
- #undef mkl_dc_getrs
- #undef mkl_dc_lapacke_getrs_convert
- #undef mkl_dc_getrs_trans
- #undef mkl_dc_getrs_nxnrhs
- #undef mkl_dc_getrs_generic
- #undef mkl_dc_geqrf
- #undef mkl_dc_lapacke_geqrf_convert
- #undef mkl_dc_geqrf_mxn
- #undef mkl_dc_geqrf_generic
- #undef mkl_dc_geqrf_hh
- #undef mkl_dc_geqrf_apply
- #undef mkl_dc_potrf
- #undef mkl_dc_lapacke_potrf_convert
- #undef mkl_dc_potrf_1
- #undef mkl_dc_potrf_2
- #undef mkl_dc_potrf_3
- #undef mkl_dc_potrf_4
- #undef mkl_dc_potrf_5
- #undef mkl_dc_potrf_n
- #undef mkl_dc_potrf_generic
- #undef mkl_dc_potrf_uplo
- #if defined(MKL_DOUBLE)
- #define mkl_dc_getrf mkl_dc_dgetrf
- #define mkl_dc_lapacke_getrf_convert mkl_dc_lapacke_dgetrf_convert
- #define mkl_dc_getrfnp mkl_dc_dgetrfnp
- #define mkl_dc_lapacke_getrfnp_convert mkl_dc_lapacke_dgetrfnp_convert
- #define mkl_dc_getri mkl_dc_dgetri
- #define mkl_dc_lapacke_getri_convert mkl_dc_lapacke_dgetri_convert
- #define mkl_dc_getrs mkl_dc_dgetrs
- #define mkl_dc_lapacke_getrs_convert mkl_dc_lapacke_dgetrs_convert
- #define mkl_dc_geqrf mkl_dc_dgeqrf
- #define mkl_dc_lapacke_geqrf_convert mkl_dc_lapacke_dgeqrf_convert
- #define mkl_dc_potrf mkl_dc_dpotrf
- #define mkl_dc_lapacke_potrf_convert mkl_dc_lapacke_dpotrf_convert
- #elif defined(MKL_SINGLE)
- #define mkl_dc_getrf mkl_dc_sgetrf
- #define mkl_dc_lapacke_getrf_convert mkl_dc_lapacke_sgetrf_convert
- #define mkl_dc_getrfnp mkl_dc_sgetrfnp
- #define mkl_dc_lapacke_getrfnp_convert mkl_dc_lapacke_sgetrfnp_convert
- #define mkl_dc_getri mkl_dc_sgetri
- #define mkl_dc_lapacke_getri_convert mkl_dc_lapacke_sgetri_convert
- #define mkl_dc_getrs mkl_dc_sgetrs
- #define mkl_dc_lapacke_getrs_convert mkl_dc_lapacke_sgetrs_convert
- #define mkl_dc_geqrf mkl_dc_sgeqrf
- #define mkl_dc_lapacke_geqrf_convert mkl_dc_lapacke_sgeqrf_convert
- #define mkl_dc_potrf mkl_dc_spotrf
- #define mkl_dc_lapacke_potrf_convert mkl_dc_lapacke_spotrf_convert
- #elif defined(MKL_COMPLEX)
- #define mkl_dc_getrf mkl_dc_cgetrf
- #define mkl_dc_lapacke_getrf_convert mkl_dc_lapacke_cgetrf_convert
- #define mkl_dc_getrfnp mkl_dc_cgetrfnp
- #define mkl_dc_lapacke_getrfnp_convert mkl_dc_lapacke_cgetrfnp_convert
- #define mkl_dc_getri mkl_dc_cgetri
- #define mkl_dc_lapacke_getri_convert mkl_dc_lapacke_cgetri_convert
- #define mkl_dc_getrs mkl_dc_cgetrs
- #define mkl_dc_lapacke_getrs_convert mkl_dc_lapacke_cgetrs_convert
- #define mkl_dc_geqrf mkl_dc_cgeqrf
- #define mkl_dc_lapacke_geqrf_convert mkl_dc_lapacke_cgeqrf_convert
- #define mkl_dc_potrf mkl_dc_cpotrf
- #define mkl_dc_lapacke_potrf_convert mkl_dc_lapacke_cpotrf_convert
- #elif defined(MKL_COMPLEX16)
- #define mkl_dc_getrf mkl_dc_zgetrf
- #define mkl_dc_lapacke_getrf_convert mkl_dc_lapacke_zgetrf_convert
- #define mkl_dc_getrfnp mkl_dc_zgetrfnp
- #define mkl_dc_lapacke_getrfnp_convert mkl_dc_lapacke_zgetrfnp_convert
- #define mkl_dc_getri mkl_dc_zgetri
- #define mkl_dc_lapacke_getri_convert mkl_dc_lapacke_zgetri_convert
- #define mkl_dc_getrs mkl_dc_zgetrs
- #define mkl_dc_lapacke_getrs_convert mkl_dc_lapacke_zgetrs_convert
- #define mkl_dc_geqrf mkl_dc_zgeqrf
- #define mkl_dc_lapacke_geqrf_convert mkl_dc_lapacke_zgeqrf_convert
- #define mkl_dc_potrf mkl_dc_zpotrf
- #define mkl_dc_lapacke_potrf_convert mkl_dc_lapacke_zpotrf_convert
- #endif
- /* ?GETRF */
- #define mkl_dc_getrf_mxn(m, n, a, lda, ipiv, info, dc_access) \
- do { \
- MKL_INT i, j, k, ii, jj; \
- *info = 0; \
- for (j = 0; j < MKL_DC_MIN(m, n); ++j) { \
- mkl_dc_real_type pivot = MKL_DC_ABS1(dc_access(a, lda, j, j)); \
- MKL_INT jp = j; \
- for (i = j + 1; i < m; ++i) { \
- mkl_dc_real_type d = MKL_DC_ABS1(dc_access(a, lda, i, j)); \
- if (d > pivot) { \
- pivot = d; \
- jp = i; \
- } \
- } \
- ipiv[j] = jp + 1; \
- if (!MKL_DC_IS_ZERO(dc_access(a, lda, jp, j))) { \
- if (jp != j) { \
- for (k = 0; k < n; k++) { \
- MKL_DC_SWAP(dc_access(a, lda, j, k), dc_access(a, lda, jp, k)); \
- } \
- } \
- for (i = j + 1; i < m; ++i) { \
- MKL_DC_DIV(dc_access(a, lda, i, j), dc_access(a, lda, i, j), dc_access(a, lda, j, j)); \
- } \
- } else if (*info == 0) { \
- *info = j + 1; \
- } \
- for (jj = j + 1; jj < n; ++jj) { \
- mkl_dc_type a_j_jj = dc_access(a, lda, j, jj); \
- MKL_DC_PRAGMA_VECTOR \
- for (ii = j + 1; ii < m; ++ii) { \
- MKL_DC_SUB_MUL(dc_access(a, lda, ii, jj), dc_access(a, lda, ii, jj), \
- dc_access(a, lda, ii, j), a_j_jj); \
- } \
- } \
- } \
- } while (0) \
- #define mkl_dc_getrf_1x1(m, n, a, lda, ipiv, info, dc_access) \
- do { \
- ipiv[0] = 1; \
- *info = !MKL_DC_IS_ZERO(dc_access(a, lda, 0, 0)) ? 0 : 1; \
- } while (0) \
- #define mkl_dc_getrf_2x2(m, n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_real_type pivot = MKL_DC_ABS1(a00); \
- mkl_dc_real_type d = MKL_DC_ABS1(a10); \
- MKL_INT jp = 0; \
- *info = 0; \
- if (d > pivot) { \
- pivot = d; \
- jp = 1; \
- } \
- ipiv[0] = jp + 1; \
- if (jp == 0) { \
- if (!MKL_DC_IS_ZERO(a00)) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- } else { \
- *info = 1; \
- } \
- } else { \
- MKL_DC_SWAP(a00, a10); \
- MKL_DC_SWAP(a01, a11); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- } \
- ipiv[1] = 2; \
- if (MKL_DC_IS_ZERO(a11) && *info == 0) { \
- *info = 2; \
- } \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 1, 1) = a11; \
- } while (0) \
- #define mkl_dc_getrf_3x3(m, n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_real_type pivot = MKL_DC_ABS1(a00); \
- mkl_dc_real_type d = MKL_DC_ABS1(a10); \
- MKL_INT jp = 0; \
- *info = 0; \
- if (d > pivot) { \
- pivot = d; \
- jp = 1; \
- } \
- d = MKL_DC_ABS1(a20); \
- if (d > pivot) { \
- pivot = d; \
- jp = 2; \
- } \
- ipiv[0] = jp + 1; \
- switch (jp) { \
- case 0: \
- if (!MKL_DC_IS_ZERO(a00)) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 1: \
- if (!MKL_DC_IS_ZERO(a10)) { \
- MKL_DC_SWAP(a00, a10); \
- MKL_DC_SWAP(a01, a11); \
- MKL_DC_SWAP(a02, a12); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 2: \
- if (!MKL_DC_IS_ZERO(a20)) { \
- MKL_DC_SWAP(a00, a20); \
- MKL_DC_SWAP(a01, a21); \
- MKL_DC_SWAP(a02, a22); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- MKL_DC_SUB_MUL(a21, a21, a20, a01); \
- MKL_DC_SUB_MUL(a12, a12, a10, a02); \
- MKL_DC_SUB_MUL(a22, a22, a20, a02); \
- pivot = MKL_DC_ABS1(a11); \
- jp = 1; \
- d = MKL_DC_ABS1(a21); \
- if (d > pivot) { \
- pivot = d; \
- jp = 2; \
- } \
- ipiv[1] = jp + 1; \
- switch (jp) { \
- case 1: \
- if (!MKL_DC_IS_ZERO(a11)) { \
- MKL_DC_DIV(a21, a21, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- case 2: \
- if (!MKL_DC_IS_ZERO(a21)) { \
- MKL_DC_SWAP(a10, a20); \
- MKL_DC_SWAP(a11, a21); \
- MKL_DC_SWAP(a12, a22); \
- MKL_DC_DIV(a21, a21, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a22, a22, a21, a12); \
- ipiv[2] = 3; \
- if (MKL_DC_IS_ZERO(a22) && *info == 0) { \
- *info = 3; \
- } \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 2, 2) = a22; \
- } while (0) \
- #define mkl_dc_getrf_4x4(m, n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- mkl_dc_real_type pivot = MKL_DC_ABS1(a00); \
- mkl_dc_real_type d = MKL_DC_ABS1(a10); \
- MKL_INT jp = 0; \
- *info = 0; \
- if (d > pivot) { \
- pivot = d; \
- jp = 1; \
- } \
- d = MKL_DC_ABS1(a20); \
- if (d > pivot) { \
- pivot = d; \
- jp = 2; \
- } \
- d = MKL_DC_ABS1(a30); \
- if (d > pivot) { \
- pivot = d; \
- jp = 3; \
- } \
- ipiv[0] = jp + 1; \
- switch (jp) { \
- case 0: \
- if (!MKL_DC_IS_ZERO(a00)) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 1: \
- if (!MKL_DC_IS_ZERO(a10)) { \
- MKL_DC_SWAP(a00, a10); \
- MKL_DC_SWAP(a01, a11); \
- MKL_DC_SWAP(a02, a12); \
- MKL_DC_SWAP(a03, a13); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 2: \
- if (!MKL_DC_IS_ZERO(a20)) { \
- MKL_DC_SWAP(a00, a20); \
- MKL_DC_SWAP(a01, a21); \
- MKL_DC_SWAP(a02, a22); \
- MKL_DC_SWAP(a03, a23); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 3: \
- if (!MKL_DC_IS_ZERO(a30)) { \
- MKL_DC_SWAP(a00, a30); \
- MKL_DC_SWAP(a01, a31); \
- MKL_DC_SWAP(a02, a32); \
- MKL_DC_SWAP(a03, a33); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- MKL_DC_SUB_MUL(a21, a21, a20, a01); \
- MKL_DC_SUB_MUL(a31, a31, a30, a01); \
- MKL_DC_SUB_MUL(a12, a12, a10, a02); \
- MKL_DC_SUB_MUL(a22, a22, a20, a02); \
- MKL_DC_SUB_MUL(a32, a32, a30, a02); \
- MKL_DC_SUB_MUL(a13, a13, a10, a03); \
- MKL_DC_SUB_MUL(a23, a23, a20, a03); \
- MKL_DC_SUB_MUL(a33, a33, a30, a03); \
- pivot = MKL_DC_ABS1(a11); \
- jp = 1; \
- d = MKL_DC_ABS1(a21); \
- if (d > pivot) { \
- pivot = d; \
- jp = 2; \
- } \
- d = MKL_DC_ABS1(a31); \
- if (d > pivot) { \
- pivot = d; \
- jp = 3; \
- } \
- ipiv[1] = jp + 1; \
- switch (jp) { \
- case 1: \
- if (!MKL_DC_IS_ZERO(a11)) { \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- case 2: \
- if (!MKL_DC_IS_ZERO(a21)) { \
- MKL_DC_SWAP(a10, a20); \
- MKL_DC_SWAP(a11, a21); \
- MKL_DC_SWAP(a12, a22); \
- MKL_DC_SWAP(a13, a23); \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- case 3: \
- if (!MKL_DC_IS_ZERO(a31)) { \
- MKL_DC_SWAP(a10, a30); \
- MKL_DC_SWAP(a11, a31); \
- MKL_DC_SWAP(a12, a32); \
- MKL_DC_SWAP(a13, a33); \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a22, a22, a21, a12); \
- MKL_DC_SUB_MUL(a32, a32, a31, a12); \
- MKL_DC_SUB_MUL(a23, a23, a21, a13); \
- MKL_DC_SUB_MUL(a33, a33, a31, a13); \
- pivot = MKL_DC_ABS1(a22); \
- jp = 2; \
- d = MKL_DC_ABS1(a32); \
- if (d > pivot) { \
- pivot = d; \
- jp = 3; \
- } \
- ipiv[2] = jp + 1; \
- switch (jp) { \
- case 2: \
- if (!MKL_DC_IS_ZERO(a22)) { \
- MKL_DC_DIV(a32, a32, a22); \
- } else if (*info == 0) { \
- *info = 3; \
- } \
- break; \
- case 3: \
- if (!MKL_DC_IS_ZERO(a32)) { \
- MKL_DC_SWAP(a20, a30); \
- MKL_DC_SWAP(a21, a31); \
- MKL_DC_SWAP(a22, a32); \
- MKL_DC_SWAP(a23, a33); \
- MKL_DC_DIV(a32, a32, a22); \
- } else if (*info == 0) { \
- *info = 3; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a33, a33, a32, a23); \
- ipiv[3] = 4; \
- if (MKL_DC_IS_ZERO(a33) && *info == 0) { \
- *info = 4; \
- } \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 2, 2) = a22; \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 0, 3) = a03; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 3, 3) = a33; \
- } while (0)
- #define mkl_dc_getrf_5x5(m, n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a40 = dc_access(a, lda, 4, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a41 = dc_access(a, lda, 4, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a42 = dc_access(a, lda, 4, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- mkl_dc_type a43 = dc_access(a, lda, 4, 3); \
- mkl_dc_type a04 = dc_access(a, lda, 0, 4); \
- mkl_dc_type a14 = dc_access(a, lda, 1, 4); \
- mkl_dc_type a24 = dc_access(a, lda, 2, 4); \
- mkl_dc_type a34 = dc_access(a, lda, 3, 4); \
- mkl_dc_type a44 = dc_access(a, lda, 4, 4); \
- mkl_dc_real_type pivot = MKL_DC_ABS1(a00); \
- mkl_dc_real_type d = MKL_DC_ABS1(a10); \
- MKL_INT jp = 0; \
- *info = 0; \
- if (d > pivot) { \
- pivot = d; \
- jp = 1; \
- } \
- d = MKL_DC_ABS1(a20); \
- if (d > pivot) { \
- pivot = d; \
- jp = 2; \
- } \
- d = MKL_DC_ABS1(a30); \
- if (d > pivot) { \
- pivot = d; \
- jp = 3; \
- } \
- d = MKL_DC_ABS1(a40); \
- if (d > pivot) { \
- pivot = d; \
- jp = 4; \
- } \
- ipiv[0] = jp + 1; \
- switch (jp) { \
- case 0: \
- if (!MKL_DC_IS_ZERO(a00)) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- MKL_DC_DIV(a40, a40, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 1: \
- if (!MKL_DC_IS_ZERO(a10)) { \
- MKL_DC_SWAP(a00, a10); \
- MKL_DC_SWAP(a01, a11); \
- MKL_DC_SWAP(a02, a12); \
- MKL_DC_SWAP(a03, a13); \
- MKL_DC_SWAP(a04, a14); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- MKL_DC_DIV(a40, a40, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 2: \
- if (!MKL_DC_IS_ZERO(a20)) { \
- MKL_DC_SWAP(a00, a20); \
- MKL_DC_SWAP(a01, a21); \
- MKL_DC_SWAP(a02, a22); \
- MKL_DC_SWAP(a03, a23); \
- MKL_DC_SWAP(a04, a24); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- MKL_DC_DIV(a40, a40, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 3: \
- if (!MKL_DC_IS_ZERO(a30)) { \
- MKL_DC_SWAP(a00, a30); \
- MKL_DC_SWAP(a01, a31); \
- MKL_DC_SWAP(a02, a32); \
- MKL_DC_SWAP(a03, a33); \
- MKL_DC_SWAP(a04, a34); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- MKL_DC_DIV(a40, a40, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- case 4: \
- if (!MKL_DC_IS_ZERO(a40)) { \
- MKL_DC_SWAP(a00, a40); \
- MKL_DC_SWAP(a01, a41); \
- MKL_DC_SWAP(a02, a42); \
- MKL_DC_SWAP(a03, a43); \
- MKL_DC_SWAP(a04, a44); \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- MKL_DC_DIV(a40, a40, a00); \
- } else { \
- *info = 1; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- MKL_DC_SUB_MUL(a21, a21, a20, a01); \
- MKL_DC_SUB_MUL(a31, a31, a30, a01); \
- MKL_DC_SUB_MUL(a41, a41, a40, a01); \
- MKL_DC_SUB_MUL(a12, a12, a10, a02); \
- MKL_DC_SUB_MUL(a22, a22, a20, a02); \
- MKL_DC_SUB_MUL(a32, a32, a30, a02); \
- MKL_DC_SUB_MUL(a42, a42, a40, a02); \
- MKL_DC_SUB_MUL(a13, a13, a10, a03); \
- MKL_DC_SUB_MUL(a23, a23, a20, a03); \
- MKL_DC_SUB_MUL(a33, a33, a30, a03); \
- MKL_DC_SUB_MUL(a43, a43, a40, a03); \
- MKL_DC_SUB_MUL(a14, a14, a10, a04); \
- MKL_DC_SUB_MUL(a24, a24, a20, a04); \
- MKL_DC_SUB_MUL(a34, a34, a30, a04); \
- MKL_DC_SUB_MUL(a44, a44, a40, a04); \
- pivot = MKL_DC_ABS1(a11); \
- jp = 1; \
- d = MKL_DC_ABS1(a21); \
- if (d > pivot) { \
- pivot = d; \
- jp = 2; \
- } \
- d = MKL_DC_ABS1(a31); \
- if (d > pivot) { \
- pivot = d; \
- jp = 3; \
- } \
- d = MKL_DC_ABS1(a41); \
- if (d > pivot) { \
- pivot = d; \
- jp = 4; \
- } \
- ipiv[1] = jp + 1; \
- switch (jp) { \
- case 1: \
- if (!MKL_DC_IS_ZERO(a11)) { \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- MKL_DC_DIV(a41, a41, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- case 2: \
- if (!MKL_DC_IS_ZERO(a21)) { \
- MKL_DC_SWAP(a10, a20); \
- MKL_DC_SWAP(a11, a21); \
- MKL_DC_SWAP(a12, a22); \
- MKL_DC_SWAP(a13, a23); \
- MKL_DC_SWAP(a14, a24); \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- MKL_DC_DIV(a41, a41, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- case 3: \
- if (!MKL_DC_IS_ZERO(a31)) { \
- MKL_DC_SWAP(a10, a30); \
- MKL_DC_SWAP(a11, a31); \
- MKL_DC_SWAP(a12, a32); \
- MKL_DC_SWAP(a13, a33); \
- MKL_DC_SWAP(a14, a34); \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- MKL_DC_DIV(a41, a41, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- case 4: \
- if (!MKL_DC_IS_ZERO(a41)) { \
- MKL_DC_SWAP(a10, a40); \
- MKL_DC_SWAP(a11, a41); \
- MKL_DC_SWAP(a12, a42); \
- MKL_DC_SWAP(a13, a43); \
- MKL_DC_SWAP(a14, a44); \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- MKL_DC_DIV(a41, a41, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a22, a22, a21, a12); \
- MKL_DC_SUB_MUL(a32, a32, a31, a12); \
- MKL_DC_SUB_MUL(a42, a42, a41, a12); \
- MKL_DC_SUB_MUL(a23, a23, a21, a13); \
- MKL_DC_SUB_MUL(a33, a33, a31, a13); \
- MKL_DC_SUB_MUL(a43, a43, a41, a13); \
- MKL_DC_SUB_MUL(a24, a24, a21, a14); \
- MKL_DC_SUB_MUL(a34, a34, a31, a14); \
- MKL_DC_SUB_MUL(a44, a44, a41, a14); \
- pivot = MKL_DC_ABS1(a22); \
- jp = 2; \
- d = MKL_DC_ABS1(a32); \
- if (d > pivot) { \
- pivot = d; \
- jp = 3; \
- } \
- d = MKL_DC_ABS1(a42); \
- if (d > pivot) { \
- pivot = d; \
- jp = 4; \
- } \
- ipiv[2] = jp + 1; \
- switch (jp) { \
- case 2: \
- if (!MKL_DC_IS_ZERO(a22)) { \
- MKL_DC_DIV(a32, a32, a22); \
- MKL_DC_DIV(a42, a42, a22); \
- } else if (*info == 0) { \
- *info = 3; \
- } \
- break; \
- case 3: \
- if (!MKL_DC_IS_ZERO(a32)) { \
- MKL_DC_SWAP(a20, a30); \
- MKL_DC_SWAP(a21, a31); \
- MKL_DC_SWAP(a22, a32); \
- MKL_DC_SWAP(a23, a33); \
- MKL_DC_SWAP(a24, a34); \
- MKL_DC_DIV(a32, a32, a22); \
- MKL_DC_DIV(a42, a42, a22); \
- } else if (*info == 0) { \
- *info = 3; \
- } \
- break; \
- case 4: \
- if (!MKL_DC_IS_ZERO(a42)) { \
- MKL_DC_SWAP(a20, a40); \
- MKL_DC_SWAP(a21, a41); \
- MKL_DC_SWAP(a22, a42); \
- MKL_DC_SWAP(a23, a43); \
- MKL_DC_SWAP(a24, a44); \
- MKL_DC_DIV(a32, a32, a22); \
- MKL_DC_DIV(a42, a42, a22); \
- } else if (*info == 0) { \
- *info = 3; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a33, a33, a32, a23); \
- MKL_DC_SUB_MUL(a43, a43, a42, a23); \
- MKL_DC_SUB_MUL(a34, a34, a32, a24); \
- MKL_DC_SUB_MUL(a44, a44, a42, a24); \
- pivot = MKL_DC_ABS1(a33); \
- jp = 3; \
- d = MKL_DC_ABS1(a43); \
- if (d > pivot) { \
- pivot = d; \
- jp = 4; \
- } \
- ipiv[3] = jp + 1; \
- switch (jp) { \
- case 3: \
- if (!MKL_DC_IS_ZERO(a33)) { \
- MKL_DC_DIV(a43, a43, a33); \
- } else if (*info == 0) { \
- *info = 4; \
- } \
- break; \
- case 4: \
- if (!MKL_DC_IS_ZERO(a43)) { \
- MKL_DC_SWAP(a30, a40); \
- MKL_DC_SWAP(a31, a41); \
- MKL_DC_SWAP(a32, a42); \
- MKL_DC_SWAP(a33, a43); \
- MKL_DC_SWAP(a34, a44); \
- MKL_DC_DIV(a43, a43, a33); \
- } else if (*info == 0) { \
- *info = 4; \
- } \
- break; \
- } \
- MKL_DC_SUB_MUL(a44, a44, a43, a34); \
- ipiv[4] = 5; \
- if (MKL_DC_IS_ZERO(a44) && *info == 0) { \
- *info = 5; \
- } \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 4, 0) = a40; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 4, 1) = a41; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 2, 2) = a22; \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 4, 2) = a42; \
- dc_access(a, lda, 0, 3) = a03; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 3, 3) = a33; \
- dc_access(a, lda, 4, 3) = a43; \
- dc_access(a, lda, 0, 4) = a04; \
- dc_access(a, lda, 1, 4) = a14; \
- dc_access(a, lda, 2, 4) = a24; \
- dc_access(a, lda, 3, 4) = a34; \
- dc_access(a, lda, 4, 4) = a44; \
- } while (0)
- #define mkl_dc_getrf_generic(m, n, a, lda, ipiv, info, dc_access) \
- do { \
- if (m == n && m <= 5) { \
- switch (m) { \
- case 0: \
- break; \
- case 1: \
- mkl_dc_getrf_1x1(m, n, a, lda, ipiv, info, dc_access); \
- break; \
- case 2: \
- mkl_dc_getrf_2x2(m, n, a, lda, ipiv, info, dc_access); \
- break; \
- case 3: \
- mkl_dc_getrf_3x3(m, n, a, lda, ipiv, info, dc_access); \
- break; \
- case 4: \
- mkl_dc_getrf_4x4(m, n, a, lda, ipiv, info, dc_access); \
- break; \
- case 5: \
- mkl_dc_getrf_5x5(m, n, a, lda, ipiv, info, dc_access); \
- break; \
- } \
- } else { \
- mkl_dc_getrf_mxn(m, n, a, lda, ipiv, info, dc_access); \
- } \
- } while (0)
- static __inline void mkl_dc_getrf(MKL_INT m, MKL_INT n, mkl_dc_type* a, MKL_INT lda, MKL_INT* ipiv, MKL_INT* info)
- {
- mkl_dc_getrf_generic(m, n, a, lda, ipiv, info, MKL_DC_MN);
- }
- #ifndef MKL_DIRECT_CALL_LAPACKE_DISABLE
- static __inline lapack_int mkl_dc_lapacke_getrf_convert(int matrix_layout, lapack_int m, lapack_int n, mkl_dc_type* a, lapack_int lda, lapack_int* ipiv)
- {
- lapack_int info = 0;
- if (MKL_DC_GETRF_CHECKSIZE(m, n)) {
- if (matrix_layout == LAPACK_ROW_MAJOR) {
- mkl_dc_getrf_generic(m, n, a, lda, ipiv, &info, MKL_DC_MT);
- } else {
- mkl_dc_getrf_generic(m, n, a, lda, ipiv, &info, MKL_DC_MN);
- }
- return info;
- }
- return MKL_DC_CONCAT3(LAPACKE_, MKL_DC_PREC_LETTER, getrf)(matrix_layout, m, n, a, lda, ipiv);
- }
- #endif
- /* ?GETRFNP */
- #define mkl_dc_getrfnp_mxn(m, n, a, lda, info, dc_access) \
- do { \
- MKL_INT i; \
- MKL_INT j; \
- MKL_INT ii; \
- MKL_INT jj; \
- MKL_INT is_zero; \
- for (j = 0; j < MKL_DC_MIN(m, n); j++) { \
- mkl_dc_type ajj = dc_access(a, lda, j, j); \
- MKL_DC_IS_ZERO2(is_zero, ajj); \
- if (!is_zero) { \
- for (i = j + 1; i < m; i++) { \
- MKL_DC_DIV(dc_access(a, lda, i, j), dc_access(a, lda, i, j), dc_access(a, lda, j, j)); \
- } \
- } else if (*info == 0) { \
- *info = j + 1; \
- } \
- for (jj = j + 1; jj < n; jj++) { \
- mkl_dc_type a_j_jj = dc_access(a, lda, j, jj); \
- for (ii = j + 1; ii < m; ii++) { \
- MKL_DC_SUB_MUL(dc_access(a, lda, ii, jj), dc_access(a, lda, ii, jj), dc_access(a, lda, ii, j), a_j_jj); \
- } \
- } \
- } \
- } while (0) \
- #define mkl_dc_getrfnp_1x1(m, n, a, lda, info, dc_access) \
- do { \
- MKL_INT is_zero; \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- MKL_DC_IS_ZERO2(is_zero, a00); \
- *info = !!is_zero; \
- } while (0) \
- #define mkl_dc_getrfnp_2x2(m, n, a, lda, info, dc_access) \
- do { \
- MKL_INT is_zero; \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- MKL_DC_IS_ZERO2(is_zero, a00); \
- if (!is_zero) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- } else { \
- *info = 1; \
- } \
- MKL_DC_IS_ZERO2(is_zero, a11); \
- if (is_zero && *info == 0) { \
- *info = 2; \
- } \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 1, 1) = a11; \
- } while (0) \
- #define mkl_dc_getrfnp_3x3(m, n, a, lda, info, dc_access) \
- do { \
- MKL_INT is_zero; \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- MKL_DC_IS_ZERO2(is_zero, a00); \
- if (!is_zero) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- } else { \
- *info = 1; \
- } \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- MKL_DC_SUB_MUL(a21, a21, a20, a01); \
- MKL_DC_SUB_MUL(a12, a12, a10, a02); \
- MKL_DC_SUB_MUL(a22, a22, a20, a02); \
- MKL_DC_IS_ZERO2(is_zero, a11); \
- if (!is_zero) { \
- MKL_DC_DIV(a21, a21, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- MKL_DC_SUB_MUL(a22, a22, a21, a12); \
- MKL_DC_IS_ZERO2(is_zero, a22); \
- if (is_zero && *info == 0) { \
- *info = 3; \
- } \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 2, 2) = a22; \
- } while (0) \
- #define mkl_dc_getrfnp_4x4(m, n, a, lda, info, dc_access) \
- do { \
- MKL_INT is_zero; \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- MKL_DC_IS_ZERO2(is_zero, a00); \
- if (!is_zero) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- } else { \
- *info = 1; \
- } \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- MKL_DC_SUB_MUL(a21, a21, a20, a01); \
- MKL_DC_SUB_MUL(a31, a31, a30, a01); \
- MKL_DC_SUB_MUL(a12, a12, a10, a02); \
- MKL_DC_SUB_MUL(a22, a22, a20, a02); \
- MKL_DC_SUB_MUL(a32, a32, a30, a02); \
- MKL_DC_SUB_MUL(a13, a13, a10, a03); \
- MKL_DC_SUB_MUL(a23, a23, a20, a03); \
- MKL_DC_SUB_MUL(a33, a33, a30, a03); \
- MKL_DC_IS_ZERO2(is_zero, a11); \
- if (!is_zero) { \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- MKL_DC_SUB_MUL(a22, a22, a21, a12); \
- MKL_DC_SUB_MUL(a32, a32, a31, a12); \
- MKL_DC_SUB_MUL(a23, a23, a21, a13); \
- MKL_DC_SUB_MUL(a33, a33, a31, a13); \
- MKL_DC_IS_ZERO2(is_zero, a22); \
- if (!is_zero) { \
- MKL_DC_DIV(a32, a32, a22); \
- } else if (*info == 0) { \
- *info = 3; \
- } \
- MKL_DC_SUB_MUL(a33, a33, a32, a23); \
- MKL_DC_IS_ZERO2(is_zero, a33); \
- if (is_zero && *info == 0) { \
- *info = 4; \
- } \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 2, 2) = a22; \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 3, 3) = a33; \
- } while (0)
- #define mkl_dc_getrfnp_5x5(m, n, a, lda, info, dc_access) \
- do { \
- MKL_INT is_zero; \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a40 = dc_access(a, lda, 4, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a41 = dc_access(a, lda, 4, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a42 = dc_access(a, lda, 4, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- mkl_dc_type a43 = dc_access(a, lda, 4, 3); \
- mkl_dc_type a04 = dc_access(a, lda, 0, 4); \
- mkl_dc_type a14 = dc_access(a, lda, 1, 4); \
- mkl_dc_type a24 = dc_access(a, lda, 2, 4); \
- mkl_dc_type a34 = dc_access(a, lda, 3, 4); \
- mkl_dc_type a44 = dc_access(a, lda, 4, 4); \
- MKL_DC_IS_ZERO2(is_zero, a00); \
- if (!is_zero) { \
- MKL_DC_DIV(a10, a10, a00); \
- MKL_DC_DIV(a20, a20, a00); \
- MKL_DC_DIV(a30, a30, a00); \
- MKL_DC_DIV(a40, a40, a00); \
- } else { \
- *info = 1; \
- } \
- MKL_DC_SUB_MUL(a11, a11, a10, a01); \
- MKL_DC_SUB_MUL(a21, a21, a20, a01); \
- MKL_DC_SUB_MUL(a31, a31, a30, a01); \
- MKL_DC_SUB_MUL(a41, a41, a40, a01); \
- MKL_DC_SUB_MUL(a12, a12, a10, a02); \
- MKL_DC_SUB_MUL(a22, a22, a20, a02); \
- MKL_DC_SUB_MUL(a32, a32, a30, a02); \
- MKL_DC_SUB_MUL(a42, a42, a40, a02); \
- MKL_DC_SUB_MUL(a13, a13, a10, a03); \
- MKL_DC_SUB_MUL(a23, a23, a20, a03); \
- MKL_DC_SUB_MUL(a33, a33, a30, a03); \
- MKL_DC_SUB_MUL(a43, a43, a40, a03); \
- MKL_DC_SUB_MUL(a14, a14, a10, a04); \
- MKL_DC_SUB_MUL(a24, a24, a20, a04); \
- MKL_DC_SUB_MUL(a34, a34, a30, a04); \
- MKL_DC_SUB_MUL(a44, a44, a40, a04); \
- MKL_DC_IS_ZERO2(is_zero, a11); \
- if (!is_zero) { \
- MKL_DC_DIV(a21, a21, a11); \
- MKL_DC_DIV(a31, a31, a11); \
- MKL_DC_DIV(a41, a41, a11); \
- } else if (*info == 0) { \
- *info = 2; \
- } \
- MKL_DC_SUB_MUL(a22, a22, a21, a12); \
- MKL_DC_SUB_MUL(a32, a32, a31, a12); \
- MKL_DC_SUB_MUL(a42, a42, a41, a12); \
- MKL_DC_SUB_MUL(a23, a23, a21, a13); \
- MKL_DC_SUB_MUL(a33, a33, a31, a13); \
- MKL_DC_SUB_MUL(a43, a43, a41, a13); \
- MKL_DC_SUB_MUL(a24, a24, a21, a14); \
- MKL_DC_SUB_MUL(a34, a34, a31, a14); \
- MKL_DC_SUB_MUL(a44, a44, a41, a14); \
- MKL_DC_IS_ZERO2(is_zero, a22); \
- if (!is_zero) { \
- MKL_DC_DIV(a32, a32, a22); \
- MKL_DC_DIV(a42, a42, a22); \
- } else if (*info == 0) { \
- *info = 3; \
- } \
- MKL_DC_SUB_MUL(a33, a33, a32, a23); \
- MKL_DC_SUB_MUL(a43, a43, a42, a23); \
- MKL_DC_SUB_MUL(a34, a34, a32, a24); \
- MKL_DC_SUB_MUL(a44, a44, a42, a24); \
- MKL_DC_IS_ZERO2(is_zero, a33); \
- if (!is_zero) { \
- MKL_DC_DIV(a43, a43, a33); \
- } else if (*info == 0) { \
- *info = 4; \
- } \
- MKL_DC_SUB_MUL(a44, a44, a43, a34); \
- MKL_DC_IS_ZERO2(is_zero, a44); \
- if (is_zero && *info == 0) { \
- *info = 5; \
- } \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 4, 0) = a40; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 4, 1) = a41; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 2, 2) = a22; \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 4, 2) = a42; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 3, 3) = a33; \
- dc_access(a, lda, 4, 3) = a43; \
- dc_access(a, lda, 1, 4) = a14; \
- dc_access(a, lda, 2, 4) = a24; \
- dc_access(a, lda, 3, 4) = a34; \
- dc_access(a, lda, 4, 4) = a44; \
- } while (0)
- #define mkl_dc_getrfnp_generic(m, n, a, lda, info, dc_access) \
- do { \
- if (m == n && m <= 5) { \
- switch (m) { \
- case 0: \
- break; \
- case 1: \
- mkl_dc_getrfnp_1x1(m, n, a, lda, info, dc_access); \
- break; \
- case 2: \
- mkl_dc_getrfnp_2x2(m, n, a, lda, info, dc_access); \
- break; \
- case 3: \
- mkl_dc_getrfnp_3x3(m, n, a, lda, info, dc_access); \
- break; \
- case 4: \
- mkl_dc_getrfnp_4x4(m, n, a, lda, info, dc_access); \
- break; \
- case 5: \
- mkl_dc_getrfnp_5x5(m, n, a, lda, info, dc_access); \
- break; \
- } \
- } else { \
- mkl_dc_getrfnp_mxn(m, n, a, lda, info, dc_access); \
- } \
- } while (0)
- static __inline void mkl_dc_getrfnp(MKL_INT m, MKL_INT n, mkl_dc_type* a, MKL_INT lda, MKL_INT* info)
- {
- *info = 0;
- mkl_dc_getrfnp_generic(m, n, a, lda, info, MKL_DC_MN);
- }
- #ifndef MKL_DIRECT_CALL_LAPACKE_DISABLE
- static __inline lapack_int mkl_dc_lapacke_getrfnp_convert(int matrix_layout, lapack_int m, lapack_int n, mkl_dc_type* a, lapack_int lda)
- {
- lapack_int info = 0;
- if (MKL_DC_GETRFNP_CHECKSIZE(m, n)) {
- if (matrix_layout == LAPACK_ROW_MAJOR) {
- mkl_dc_getrfnp_generic(m, n, a, lda, &info, MKL_DC_MT);
- } else {
- mkl_dc_getrfnp_generic(m, n, a, lda, &info, MKL_DC_MN);
- }
- return info;
- }
- return MKL_DC_CONCAT3(LAPACKE_mkl_, MKL_DC_PREC_LETTER, getrfnp)(matrix_layout, m, n, a, lda);
- }
- #endif
- /* ?GETRS */
- #define mkl_dc_getrs_nxnrhs(trans, n, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access) \
- do { \
- MKL_INT i, j, k; \
- if (trans == 'N') { \
- for (i = 0; i < n; ++i) { \
- MKL_INT ip = ipiv[i] - 1; \
- if (ip != i) { \
- for (k = 0; k < nrhs; k++) { \
- MKL_DC_SWAP(b_access(b, ldb, i, k), b_access(b, ldb, ip, k)); \
- } \
- } \
- } \
- for (j = 0; j < nrhs; ++j) { \
- for (k = 0; k < n; k++) { \
- if (!MKL_DC_IS_ZERO(b_access(b, ldb, k, j))) { \
- for (i = k + 1; i < n; ++i) { \
- MKL_DC_SUB_MUL(b_access(b, ldb, i, j), b_access(b, ldb, i, j), b_access(b, ldb, k, j), dc_access(a, lda, i, k)); \
- } \
- } \
- } \
- } \
- for (j = 0; j < nrhs; ++j) { \
- for (k = n - 1; k >= 0; k--) { \
- if (!MKL_DC_IS_ZERO(b_access(b, ldb, k, j))) { \
- MKL_DC_DIV(b_access(b, ldb, k, j), b_access(b, ldb, k, j), dc_access(a, lda, k, k)); \
- for (i = 0; i < k; ++i) { \
- MKL_DC_SUB_MUL(b_access(b, ldb, i, j), b_access(b, ldb, i, j), b_access(b, ldb, k, j), dc_access(a, lda, i, k)); \
- } \
- } \
- } \
- } \
- } else { \
- for (j = 0; j < nrhs; ++j) { \
- for (i = 0; i < n; ++i) { \
- mkl_dc_type temp = b_access(b, ldb, i, j); \
- mkl_dc_type aii = dc_access(a, lda, i, i); \
- for (k = 0; k < i; k++) { \
- mkl_dc_type aki = dc_access(a, lda, k, i); \
- if (trans == 'C') { \
- MKL_DC_CONJ(aki, aki); \
- } \
- MKL_DC_SUB_MUL(temp, temp, aki, b_access(b, ldb, k, j)); \
- } \
- if (trans == 'C') { \
- MKL_DC_CONJ(aii, aii); \
- } \
- MKL_DC_DIV(b_access(b, ldb, i, j), temp, aii); \
- } \
- } \
- for (j = 0; j < nrhs; ++j) { \
- for (i = n - 1; i >= 0; i--) { \
- mkl_dc_type temp = b_access(b, ldb, i, j); \
- for (k = i + 1; k < n; k++) { \
- mkl_dc_type aki = dc_access(a, lda, k, i); \
- if (trans == 'C') { \
- MKL_DC_CONJ(aki, aki); \
- } \
- MKL_DC_SUB_MUL(temp, temp, aki, b_access(b, ldb, k, j)); \
- } \
- b_access(b, ldb, i, j) = temp; \
- } \
- } \
- for (i = n - 1; i >= 0; i--) { \
- MKL_INT ip = ipiv[i] - 1; \
- if (ip != i) { \
- for (k = 0; k < nrhs; k++) { \
- MKL_DC_SWAP(b_access(b, ldb, i, k), b_access(b, ldb, ip, k)); \
- } \
- } \
- } \
- } \
- } while (0)
- #define mkl_dc_getrs_trans(trans, n, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access) \
- do { \
- switch (n) { \
- case 1: \
- mkl_dc_getrs_nxnrhs(trans, 1, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- break; \
- case 2: \
- mkl_dc_getrs_nxnrhs(trans, 2, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- break; \
- case 3: \
- mkl_dc_getrs_nxnrhs(trans, 3, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- break; \
- case 4: \
- mkl_dc_getrs_nxnrhs(trans, 4, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- break; \
- case 5: \
- mkl_dc_getrs_nxnrhs(trans, 5, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- break; \
- default: \
- mkl_dc_getrs_nxnrhs(trans, n, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- break; \
- } \
- } while (0)
- #define mkl_dc_getrs_generic(trans, n, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access) \
- do { \
- int ntrans = MKL_DC_MisN(trans); \
- int ttrans = MKL_DC_MisT(trans); \
- *info = 0; \
- if (n != 0 && nrhs != 0) { \
- if (ntrans) { \
- mkl_dc_getrs_trans('N', n, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- } else if (ttrans) { \
- mkl_dc_getrs_trans('T', n, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- } else { \
- mkl_dc_getrs_trans('C', n, nrhs, a, lda, ipiv, b, ldb, info, dc_access, b_access); \
- } \
- } \
- } while (0)
- static __inline void mkl_dc_getrs(char trans, MKL_INT n, MKL_INT nrhs, const mkl_dc_type* a, MKL_INT lda,
- const MKL_INT* ipiv, mkl_dc_type* b, MKL_INT ldb, MKL_INT* info)
- {
- mkl_dc_getrs_generic(trans, n, nrhs, a, lda, ipiv, b, ldb, info, MKL_DC_MN, MKL_DC_MN);
- }
- #ifndef MKL_DIRECT_CALL_LAPACKE_DISABLE
- static __inline lapack_int mkl_dc_lapacke_getrs_convert(int matrix_layout, char trans, lapack_int n, lapack_int nrhs,
- const mkl_dc_type* a, lapack_int lda, const lapack_int* ipiv,
- mkl_dc_type* b, lapack_int ldb)
- {
- lapack_int info = 0;
- if (MKL_DC_GETRS_CHECKSIZE(n, nrhs)) {
- if (matrix_layout == LAPACK_ROW_MAJOR) {
- mkl_dc_getrs_generic(trans, n, nrhs, a, lda, ipiv, b, ldb, &info, MKL_DC_MT, MKL_DC_MT);
- } else {
- mkl_dc_getrs_generic(trans, n, nrhs, a, lda, ipiv, b, ldb, &info, MKL_DC_MN, MKL_DC_MN);
- }
- return info;
- }
- return MKL_DC_CONCAT3(LAPACKE_, MKL_DC_PREC_LETTER, getrs)(matrix_layout, trans, n, nrhs, a, lda, ipiv, b, ldb);
- }
- #endif
- /* ?GETRI */
- #define mkl_dc_getri_1x1(n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type ajj; \
- if (MKL_DC_IS_ZERO(a00)) { \
- *info = 1; \
- break; \
- } \
- *info = 0; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_NEG(ajj, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- } while (0)
- #define mkl_dc_getri_2x2(n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type w1; \
- mkl_dc_type tmp; \
- MKL_INT jp; \
- mkl_dc_type ajj; \
- if (MKL_DC_IS_ZERO(a00)) { \
- *info = 1; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a11)) { \
- *info = 2; \
- break; \
- } \
- *info = 0; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_NEG(ajj, a00); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_NEG(ajj, a11); \
- tmp = a01; \
- MKL_DC_MUL(a01, a01, a00); \
- MKL_DC_MUL(a01, a01, ajj); \
- w1 = a10; \
- MKL_DC_SET_ZERO(a10); \
- MKL_DC_NEG(tmp, w1); \
- MKL_DC_MUL_ADD(a00, tmp, a01, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a11, a10); \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 1, 1) = a11; \
- jp = ipiv[0] - 1; \
- if (jp != 0) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 0), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 0), dc_access(a, lda, 1, jp)); \
- } \
- } while (0)
- #define mkl_dc_getri_3x3(n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type w1; \
- mkl_dc_type w2; \
- mkl_dc_type tmp; \
- MKL_INT jp; \
- mkl_dc_type ajj; \
- if (MKL_DC_IS_ZERO(a00)) { \
- *info = 1; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a11)) { \
- *info = 2; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a22)) { \
- *info = 3; \
- break; \
- } \
- *info = 0; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_NEG(ajj, a00); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_NEG(ajj, a11); \
- tmp = a01; \
- MKL_DC_MUL(a01, a01, a00); \
- MKL_DC_MUL(a01, a01, ajj); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_NEG(ajj, a22); \
- tmp = a02; \
- MKL_DC_MUL(a02, a02, a00); \
- tmp = a12; \
- MKL_DC_MUL_ADD(a02, tmp, a01, a02); \
- MKL_DC_MUL(a12, a12, a11); \
- MKL_DC_MUL(a02, a02, ajj); \
- MKL_DC_MUL(a12, a12, ajj); \
- w2 = a21; \
- MKL_DC_SET_ZERO(a21); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a01, tmp, a02, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a12, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a22, a21); \
- w1 = a10; \
- MKL_DC_SET_ZERO(a10); \
- w2 = a20; \
- MKL_DC_SET_ZERO(a20); \
- MKL_DC_NEG(tmp, w1); \
- MKL_DC_MUL_ADD(a00, tmp, a01, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a11, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a21, a20); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a00, tmp, a02, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a12, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a22, a20); \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 2, 2) = a22; \
- jp = ipiv[1] - 1; \
- if (jp != 1) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 1), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 1), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 1), dc_access(a, lda, 2, jp)); \
- } \
- jp = ipiv[0] - 1; \
- if (jp != 0) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 0), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 0), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 0), dc_access(a, lda, 2, jp)); \
- } \
- } while (0)
- #define mkl_dc_getri_4x4(n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- mkl_dc_type w1; \
- mkl_dc_type w2; \
- mkl_dc_type w3; \
- mkl_dc_type tmp; \
- MKL_INT jp; \
- mkl_dc_type ajj; \
- if (MKL_DC_IS_ZERO(a00)) { \
- *info = 1; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a11)) { \
- *info = 2; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a22)) { \
- *info = 3; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a33)) { \
- *info = 4; \
- break; \
- } \
- *info = 0; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_NEG(ajj, a00); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_NEG(ajj, a11); \
- tmp = a01; \
- MKL_DC_MUL(a01, a01, a00); \
- MKL_DC_MUL(a01, a01, ajj); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_NEG(ajj, a22); \
- tmp = a02; \
- MKL_DC_MUL(a02, a02, a00); \
- tmp = a12; \
- MKL_DC_MUL_ADD(a02, tmp, a01, a02); \
- MKL_DC_MUL(a12, a12, a11); \
- MKL_DC_MUL(a02, a02, ajj); \
- MKL_DC_MUL(a12, a12, ajj); \
- MKL_DC_INV(a33, a33); \
- MKL_DC_NEG(ajj, a33); \
- tmp = a03; \
- MKL_DC_MUL(a03, a03, a00); \
- tmp = a13; \
- MKL_DC_MUL_ADD(a03, tmp, a01, a03); \
- MKL_DC_MUL(a13, a13, a11); \
- tmp = a23; \
- MKL_DC_MUL_ADD(a03, tmp, a02, a03); \
- MKL_DC_MUL_ADD(a13, tmp, a12, a13); \
- MKL_DC_MUL(a23, a23, a22); \
- MKL_DC_MUL(a03, a03, ajj); \
- MKL_DC_MUL(a13, a13, ajj); \
- MKL_DC_MUL(a23, a23, ajj); \
- w3 = a32; \
- MKL_DC_SET_ZERO(a32); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a02, tmp, a03, a02); \
- MKL_DC_MUL_ADD(a12, tmp, a13, a12); \
- MKL_DC_MUL_ADD(a22, tmp, a23, a22); \
- MKL_DC_MUL_ADD(a32, tmp, a33, a32); \
- w2 = a21; \
- MKL_DC_SET_ZERO(a21); \
- w3 = a31; \
- MKL_DC_SET_ZERO(a31); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a01, tmp, a02, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a12, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a22, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a32, a31); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a01, tmp, a03, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a13, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a23, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a33, a31); \
- w1 = a10; \
- MKL_DC_SET_ZERO(a10); \
- w2 = a20; \
- MKL_DC_SET_ZERO(a20); \
- w3 = a30; \
- MKL_DC_SET_ZERO(a30); \
- MKL_DC_NEG(tmp, w1); \
- MKL_DC_MUL_ADD(a00, tmp, a01, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a11, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a21, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a31, a30); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a00, tmp, a02, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a12, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a22, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a32, a30); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a00, tmp, a03, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a13, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a23, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a33, a30); \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 0, 3) = a03; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 2, 2) = a22; \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 3, 3) = a33; \
- jp = ipiv[2] - 1; \
- if (jp != 2) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 2), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 2), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 2), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 2), dc_access(a, lda, 3, jp)); \
- } \
- jp = ipiv[1] - 1; \
- if (jp != 1) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 1), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 1), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 1), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 1), dc_access(a, lda, 3, jp)); \
- } \
- jp = ipiv[0] - 1; \
- if (jp != 0) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 0), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 0), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 0), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 0), dc_access(a, lda, 3, jp)); \
- } \
- } while (0)
- #define mkl_dc_getri_5x5(n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a04 = dc_access(a, lda, 0, 4); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a14 = dc_access(a, lda, 1, 4); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a24 = dc_access(a, lda, 2, 4); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- mkl_dc_type a34 = dc_access(a, lda, 3, 4); \
- mkl_dc_type a40 = dc_access(a, lda, 4, 0); \
- mkl_dc_type a41 = dc_access(a, lda, 4, 1); \
- mkl_dc_type a42 = dc_access(a, lda, 4, 2); \
- mkl_dc_type a43 = dc_access(a, lda, 4, 3); \
- mkl_dc_type a44 = dc_access(a, lda, 4, 4); \
- mkl_dc_type w1; \
- mkl_dc_type w2; \
- mkl_dc_type w3; \
- mkl_dc_type w4; \
- mkl_dc_type tmp; \
- MKL_INT jp; \
- mkl_dc_type ajj; \
- if (MKL_DC_IS_ZERO(a00)) { \
- *info = 1; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a11)) { \
- *info = 2; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a22)) { \
- *info = 3; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a33)) { \
- *info = 4; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a44)) { \
- *info = 5; \
- break; \
- } \
- *info = 0; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_NEG(ajj, a00); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_NEG(ajj, a11); \
- tmp = a01; \
- MKL_DC_MUL(a01, a01, a00); \
- MKL_DC_MUL(a01, a01, ajj); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_NEG(ajj, a22); \
- tmp = a02; \
- MKL_DC_MUL(a02, a02, a00); \
- tmp = a12; \
- MKL_DC_MUL_ADD(a02, tmp, a01, a02); \
- MKL_DC_MUL(a12, a12, a11); \
- MKL_DC_MUL(a02, a02, ajj); \
- MKL_DC_MUL(a12, a12, ajj); \
- MKL_DC_INV(a33, a33); \
- MKL_DC_NEG(ajj, a33); \
- tmp = a03; \
- MKL_DC_MUL(a03, a03, a00); \
- tmp = a13; \
- MKL_DC_MUL_ADD(a03, tmp, a01, a03); \
- MKL_DC_MUL(a13, a13, a11); \
- tmp = a23; \
- MKL_DC_MUL_ADD(a03, tmp, a02, a03); \
- MKL_DC_MUL_ADD(a13, tmp, a12, a13); \
- MKL_DC_MUL(a23, a23, a22); \
- MKL_DC_MUL(a03, a03, ajj); \
- MKL_DC_MUL(a13, a13, ajj); \
- MKL_DC_MUL(a23, a23, ajj); \
- MKL_DC_INV(a44, a44); \
- MKL_DC_NEG(ajj, a44); \
- tmp = a04; \
- MKL_DC_MUL(a04, a04, a00); \
- tmp = a14; \
- MKL_DC_MUL_ADD(a04, tmp, a01, a04); \
- MKL_DC_MUL(a14, a14, a11); \
- tmp = a24; \
- MKL_DC_MUL_ADD(a04, tmp, a02, a04); \
- MKL_DC_MUL_ADD(a14, tmp, a12, a14); \
- MKL_DC_MUL(a24, a24, a22); \
- tmp = a34; \
- MKL_DC_MUL_ADD(a04, tmp, a03, a04); \
- MKL_DC_MUL_ADD(a14, tmp, a13, a14); \
- MKL_DC_MUL_ADD(a24, tmp, a23, a24); \
- MKL_DC_MUL(a34, a34, a33); \
- MKL_DC_MUL(a04, a04, ajj); \
- MKL_DC_MUL(a14, a14, ajj); \
- MKL_DC_MUL(a24, a24, ajj); \
- MKL_DC_MUL(a34, a34, ajj); \
- w4 = a43; \
- MKL_DC_SET_ZERO(a43); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a03, tmp, a04, a03); \
- MKL_DC_MUL_ADD(a13, tmp, a14, a13); \
- MKL_DC_MUL_ADD(a23, tmp, a24, a23); \
- MKL_DC_MUL_ADD(a33, tmp, a34, a33); \
- MKL_DC_MUL_ADD(a43, tmp, a44, a43); \
- w3 = a32; \
- MKL_DC_SET_ZERO(a32); \
- w4 = a42; \
- MKL_DC_SET_ZERO(a42); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a02, tmp, a03, a02); \
- MKL_DC_MUL_ADD(a12, tmp, a13, a12); \
- MKL_DC_MUL_ADD(a22, tmp, a23, a22); \
- MKL_DC_MUL_ADD(a32, tmp, a33, a32); \
- MKL_DC_MUL_ADD(a42, tmp, a43, a42); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a02, tmp, a04, a02); \
- MKL_DC_MUL_ADD(a12, tmp, a14, a12); \
- MKL_DC_MUL_ADD(a22, tmp, a24, a22); \
- MKL_DC_MUL_ADD(a32, tmp, a34, a32); \
- MKL_DC_MUL_ADD(a42, tmp, a44, a42); \
- w2 = a21; \
- MKL_DC_SET_ZERO(a21); \
- w3 = a31; \
- MKL_DC_SET_ZERO(a31); \
- w4 = a41; \
- MKL_DC_SET_ZERO(a41); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a01, tmp, a02, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a12, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a22, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a32, a31); \
- MKL_DC_MUL_ADD(a41, tmp, a42, a41); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a01, tmp, a03, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a13, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a23, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a33, a31); \
- MKL_DC_MUL_ADD(a41, tmp, a43, a41); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a01, tmp, a04, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a14, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a24, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a34, a31); \
- MKL_DC_MUL_ADD(a41, tmp, a44, a41); \
- w1 = a10; \
- MKL_DC_SET_ZERO(a10); \
- w2 = a20; \
- MKL_DC_SET_ZERO(a20); \
- w3 = a30; \
- MKL_DC_SET_ZERO(a30); \
- w4 = a40; \
- MKL_DC_SET_ZERO(a40); \
- MKL_DC_NEG(tmp, w1); \
- MKL_DC_MUL_ADD(a00, tmp, a01, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a11, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a21, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a31, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a41, a40); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a00, tmp, a02, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a12, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a22, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a32, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a42, a40); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a00, tmp, a03, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a13, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a23, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a33, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a43, a40); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a00, tmp, a04, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a14, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a24, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a34, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a44, a40); \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 0, 3) = a03; \
- dc_access(a, lda, 0, 4) = a04; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 1, 4) = a14; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 2, 2) = a22; \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 2, 4) = a24; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 3, 3) = a33; \
- dc_access(a, lda, 3, 4) = a34; \
- dc_access(a, lda, 4, 0) = a40; \
- dc_access(a, lda, 4, 1) = a41; \
- dc_access(a, lda, 4, 2) = a42; \
- dc_access(a, lda, 4, 3) = a43; \
- dc_access(a, lda, 4, 4) = a44; \
- jp = ipiv[3] - 1; \
- if (jp != 3) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 3), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 3), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 3), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 3), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 3), dc_access(a, lda, 4, jp)); \
- } \
- jp = ipiv[2] - 1; \
- if (jp != 2) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 2), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 2), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 2), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 2), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 2), dc_access(a, lda, 4, jp)); \
- } \
- jp = ipiv[1] - 1; \
- if (jp != 1) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 1), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 1), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 1), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 1), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 1), dc_access(a, lda, 4, jp)); \
- } \
- jp = ipiv[0] - 1; \
- if (jp != 0) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 0), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 0), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 0), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 0), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 0), dc_access(a, lda, 4, jp)); \
- } \
- } while (0)
- #define mkl_dc_getri_6x6(n, a, lda, ipiv, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a04 = dc_access(a, lda, 0, 4); \
- mkl_dc_type a05 = dc_access(a, lda, 0, 5); \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a14 = dc_access(a, lda, 1, 4); \
- mkl_dc_type a15 = dc_access(a, lda, 1, 5); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a24 = dc_access(a, lda, 2, 4); \
- mkl_dc_type a25 = dc_access(a, lda, 2, 5); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- mkl_dc_type a34 = dc_access(a, lda, 3, 4); \
- mkl_dc_type a35 = dc_access(a, lda, 3, 5); \
- mkl_dc_type a40 = dc_access(a, lda, 4, 0); \
- mkl_dc_type a41 = dc_access(a, lda, 4, 1); \
- mkl_dc_type a42 = dc_access(a, lda, 4, 2); \
- mkl_dc_type a43 = dc_access(a, lda, 4, 3); \
- mkl_dc_type a44 = dc_access(a, lda, 4, 4); \
- mkl_dc_type a45 = dc_access(a, lda, 4, 5); \
- mkl_dc_type a50 = dc_access(a, lda, 5, 0); \
- mkl_dc_type a51 = dc_access(a, lda, 5, 1); \
- mkl_dc_type a52 = dc_access(a, lda, 5, 2); \
- mkl_dc_type a53 = dc_access(a, lda, 5, 3); \
- mkl_dc_type a54 = dc_access(a, lda, 5, 4); \
- mkl_dc_type a55 = dc_access(a, lda, 5, 5); \
- mkl_dc_type w1; \
- mkl_dc_type w2; \
- mkl_dc_type w3; \
- mkl_dc_type w4; \
- mkl_dc_type w5; \
- mkl_dc_type tmp; \
- MKL_INT jp; \
- mkl_dc_type ajj; \
- if (MKL_DC_IS_ZERO(a00)) { \
- *info = 1; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a11)) { \
- *info = 2; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a22)) { \
- *info = 3; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a33)) { \
- *info = 4; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a44)) { \
- *info = 5; \
- break; \
- } \
- if (MKL_DC_IS_ZERO(a55)) { \
- *info = 6; \
- break; \
- } \
- *info = 0; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_NEG(ajj, a00); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_NEG(ajj, a11); \
- tmp = a01; \
- MKL_DC_MUL(a01, a01, a00); \
- MKL_DC_MUL(a01, a01, ajj); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_NEG(ajj, a22); \
- tmp = a02; \
- MKL_DC_MUL(a02, a02, a00); \
- tmp = a12; \
- MKL_DC_MUL_ADD(a02, tmp, a01, a02); \
- MKL_DC_MUL(a12, a12, a11); \
- MKL_DC_MUL(a02, a02, ajj); \
- MKL_DC_MUL(a12, a12, ajj); \
- MKL_DC_INV(a33, a33); \
- MKL_DC_NEG(ajj, a33); \
- tmp = a03; \
- MKL_DC_MUL(a03, a03, a00); \
- tmp = a13; \
- MKL_DC_MUL_ADD(a03, tmp, a01, a03); \
- MKL_DC_MUL(a13, a13, a11); \
- tmp = a23; \
- MKL_DC_MUL_ADD(a03, tmp, a02, a03); \
- MKL_DC_MUL_ADD(a13, tmp, a12, a13); \
- MKL_DC_MUL(a23, a23, a22); \
- MKL_DC_MUL(a03, a03, ajj); \
- MKL_DC_MUL(a13, a13, ajj); \
- MKL_DC_MUL(a23, a23, ajj); \
- MKL_DC_INV(a44, a44); \
- MKL_DC_NEG(ajj, a44); \
- tmp = a04; \
- MKL_DC_MUL(a04, a04, a00); \
- tmp = a14; \
- MKL_DC_MUL_ADD(a04, tmp, a01, a04); \
- MKL_DC_MUL(a14, a14, a11); \
- tmp = a24; \
- MKL_DC_MUL_ADD(a04, tmp, a02, a04); \
- MKL_DC_MUL_ADD(a14, tmp, a12, a14); \
- MKL_DC_MUL(a24, a24, a22); \
- tmp = a34; \
- MKL_DC_MUL_ADD(a04, tmp, a03, a04); \
- MKL_DC_MUL_ADD(a14, tmp, a13, a14); \
- MKL_DC_MUL_ADD(a24, tmp, a23, a24); \
- MKL_DC_MUL(a34, a34, a33); \
- MKL_DC_MUL(a04, a04, ajj); \
- MKL_DC_MUL(a14, a14, ajj); \
- MKL_DC_MUL(a24, a24, ajj); \
- MKL_DC_MUL(a34, a34, ajj); \
- MKL_DC_INV(a55, a55); \
- MKL_DC_NEG(ajj, a55); \
- tmp = a05; \
- MKL_DC_MUL(a05, a05, a00); \
- tmp = a15; \
- MKL_DC_MUL_ADD(a05, tmp, a01, a05); \
- MKL_DC_MUL(a15, a15, a11); \
- tmp = a25; \
- MKL_DC_MUL_ADD(a05, tmp, a02, a05); \
- MKL_DC_MUL_ADD(a15, tmp, a12, a15); \
- MKL_DC_MUL(a25, a25, a22); \
- tmp = a35; \
- MKL_DC_MUL_ADD(a05, tmp, a03, a05); \
- MKL_DC_MUL_ADD(a15, tmp, a13, a15); \
- MKL_DC_MUL_ADD(a25, tmp, a23, a25); \
- MKL_DC_MUL(a35, a35, a33); \
- tmp = a45; \
- MKL_DC_MUL_ADD(a05, tmp, a04, a05); \
- MKL_DC_MUL_ADD(a15, tmp, a14, a15); \
- MKL_DC_MUL_ADD(a25, tmp, a24, a25); \
- MKL_DC_MUL_ADD(a35, tmp, a34, a35); \
- MKL_DC_MUL(a45, a45, a44); \
- MKL_DC_MUL(a05, a05, ajj); \
- MKL_DC_MUL(a15, a15, ajj); \
- MKL_DC_MUL(a25, a25, ajj); \
- MKL_DC_MUL(a35, a35, ajj); \
- MKL_DC_MUL(a45, a45, ajj); \
- w5 = a54; \
- MKL_DC_SET_ZERO(a54); \
- MKL_DC_NEG(tmp, w5); \
- MKL_DC_MUL_ADD(a04, tmp, a05, a04); \
- MKL_DC_MUL_ADD(a14, tmp, a15, a14); \
- MKL_DC_MUL_ADD(a24, tmp, a25, a24); \
- MKL_DC_MUL_ADD(a34, tmp, a35, a34); \
- MKL_DC_MUL_ADD(a44, tmp, a45, a44); \
- MKL_DC_MUL_ADD(a54, tmp, a55, a54); \
- w4 = a43; \
- MKL_DC_SET_ZERO(a43); \
- w5 = a53; \
- MKL_DC_SET_ZERO(a53); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a03, tmp, a04, a03); \
- MKL_DC_MUL_ADD(a13, tmp, a14, a13); \
- MKL_DC_MUL_ADD(a23, tmp, a24, a23); \
- MKL_DC_MUL_ADD(a33, tmp, a34, a33); \
- MKL_DC_MUL_ADD(a43, tmp, a44, a43); \
- MKL_DC_MUL_ADD(a53, tmp, a54, a53); \
- MKL_DC_NEG(tmp, w5); \
- MKL_DC_MUL_ADD(a03, tmp, a05, a03); \
- MKL_DC_MUL_ADD(a13, tmp, a15, a13); \
- MKL_DC_MUL_ADD(a23, tmp, a25, a23); \
- MKL_DC_MUL_ADD(a33, tmp, a35, a33); \
- MKL_DC_MUL_ADD(a43, tmp, a45, a43); \
- MKL_DC_MUL_ADD(a53, tmp, a55, a53); \
- w3 = a32; \
- MKL_DC_SET_ZERO(a32); \
- w4 = a42; \
- MKL_DC_SET_ZERO(a42); \
- w5 = a52; \
- MKL_DC_SET_ZERO(a52); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a02, tmp, a03, a02); \
- MKL_DC_MUL_ADD(a12, tmp, a13, a12); \
- MKL_DC_MUL_ADD(a22, tmp, a23, a22); \
- MKL_DC_MUL_ADD(a32, tmp, a33, a32); \
- MKL_DC_MUL_ADD(a42, tmp, a43, a42); \
- MKL_DC_MUL_ADD(a52, tmp, a53, a52); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a02, tmp, a04, a02); \
- MKL_DC_MUL_ADD(a12, tmp, a14, a12); \
- MKL_DC_MUL_ADD(a22, tmp, a24, a22); \
- MKL_DC_MUL_ADD(a32, tmp, a34, a32); \
- MKL_DC_MUL_ADD(a42, tmp, a44, a42); \
- MKL_DC_MUL_ADD(a52, tmp, a54, a52); \
- MKL_DC_NEG(tmp, w5); \
- MKL_DC_MUL_ADD(a02, tmp, a05, a02); \
- MKL_DC_MUL_ADD(a12, tmp, a15, a12); \
- MKL_DC_MUL_ADD(a22, tmp, a25, a22); \
- MKL_DC_MUL_ADD(a32, tmp, a35, a32); \
- MKL_DC_MUL_ADD(a42, tmp, a45, a42); \
- MKL_DC_MUL_ADD(a52, tmp, a55, a52); \
- w2 = a21; \
- MKL_DC_SET_ZERO(a21); \
- w3 = a31; \
- MKL_DC_SET_ZERO(a31); \
- w4 = a41; \
- MKL_DC_SET_ZERO(a41); \
- w5 = a51; \
- MKL_DC_SET_ZERO(a51); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a01, tmp, a02, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a12, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a22, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a32, a31); \
- MKL_DC_MUL_ADD(a41, tmp, a42, a41); \
- MKL_DC_MUL_ADD(a51, tmp, a52, a51); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a01, tmp, a03, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a13, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a23, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a33, a31); \
- MKL_DC_MUL_ADD(a41, tmp, a43, a41); \
- MKL_DC_MUL_ADD(a51, tmp, a53, a51); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a01, tmp, a04, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a14, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a24, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a34, a31); \
- MKL_DC_MUL_ADD(a41, tmp, a44, a41); \
- MKL_DC_MUL_ADD(a51, tmp, a54, a51); \
- MKL_DC_NEG(tmp, w5); \
- MKL_DC_MUL_ADD(a01, tmp, a05, a01); \
- MKL_DC_MUL_ADD(a11, tmp, a15, a11); \
- MKL_DC_MUL_ADD(a21, tmp, a25, a21); \
- MKL_DC_MUL_ADD(a31, tmp, a35, a31); \
- MKL_DC_MUL_ADD(a41, tmp, a45, a41); \
- MKL_DC_MUL_ADD(a51, tmp, a55, a51); \
- w1 = a10; \
- MKL_DC_SET_ZERO(a10); \
- w2 = a20; \
- MKL_DC_SET_ZERO(a20); \
- w3 = a30; \
- MKL_DC_SET_ZERO(a30); \
- w4 = a40; \
- MKL_DC_SET_ZERO(a40); \
- w5 = a50; \
- MKL_DC_SET_ZERO(a50); \
- MKL_DC_NEG(tmp, w1); \
- MKL_DC_MUL_ADD(a00, tmp, a01, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a11, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a21, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a31, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a41, a40); \
- MKL_DC_MUL_ADD(a50, tmp, a51, a50); \
- MKL_DC_NEG(tmp, w2); \
- MKL_DC_MUL_ADD(a00, tmp, a02, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a12, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a22, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a32, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a42, a40); \
- MKL_DC_MUL_ADD(a50, tmp, a52, a50); \
- MKL_DC_NEG(tmp, w3); \
- MKL_DC_MUL_ADD(a00, tmp, a03, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a13, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a23, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a33, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a43, a40); \
- MKL_DC_MUL_ADD(a50, tmp, a53, a50); \
- MKL_DC_NEG(tmp, w4); \
- MKL_DC_MUL_ADD(a00, tmp, a04, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a14, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a24, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a34, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a44, a40); \
- MKL_DC_MUL_ADD(a50, tmp, a54, a50); \
- MKL_DC_NEG(tmp, w5); \
- MKL_DC_MUL_ADD(a00, tmp, a05, a00); \
- MKL_DC_MUL_ADD(a10, tmp, a15, a10); \
- MKL_DC_MUL_ADD(a20, tmp, a25, a20); \
- MKL_DC_MUL_ADD(a30, tmp, a35, a30); \
- MKL_DC_MUL_ADD(a40, tmp, a45, a40); \
- MKL_DC_MUL_ADD(a50, tmp, a55, a50); \
- dc_access(a, lda, 0, 0) = a00; \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 0, 3) = a03; \
- dc_access(a, lda, 0, 4) = a04; \
- dc_access(a, lda, 0, 5) = a05; \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 1, 1) = a11; \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 1, 4) = a14; \
- dc_access(a, lda, 1, 5) = a15; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 2, 2) = a22; \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 2, 4) = a24; \
- dc_access(a, lda, 2, 5) = a25; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 3, 3) = a33; \
- dc_access(a, lda, 3, 4) = a34; \
- dc_access(a, lda, 3, 5) = a35; \
- dc_access(a, lda, 4, 0) = a40; \
- dc_access(a, lda, 4, 1) = a41; \
- dc_access(a, lda, 4, 2) = a42; \
- dc_access(a, lda, 4, 3) = a43; \
- dc_access(a, lda, 4, 4) = a44; \
- dc_access(a, lda, 4, 5) = a45; \
- dc_access(a, lda, 5, 0) = a50; \
- dc_access(a, lda, 5, 1) = a51; \
- dc_access(a, lda, 5, 2) = a52; \
- dc_access(a, lda, 5, 3) = a53; \
- dc_access(a, lda, 5, 4) = a54; \
- dc_access(a, lda, 5, 5) = a55; \
- jp = ipiv[4] - 1; \
- if (jp != 4) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 4), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 4), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 4), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 4), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 4), dc_access(a, lda, 4, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 5, 4), dc_access(a, lda, 5, jp)); \
- } \
- jp = ipiv[3] - 1; \
- if (jp != 3) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 3), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 3), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 3), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 3), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 3), dc_access(a, lda, 4, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 5, 3), dc_access(a, lda, 5, jp)); \
- } \
- jp = ipiv[2] - 1; \
- if (jp != 2) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 2), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 2), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 2), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 2), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 2), dc_access(a, lda, 4, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 5, 2), dc_access(a, lda, 5, jp)); \
- } \
- jp = ipiv[1] - 1; \
- if (jp != 1) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 1), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 1), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 1), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 1), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 1), dc_access(a, lda, 4, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 5, 1), dc_access(a, lda, 5, jp)); \
- } \
- jp = ipiv[0] - 1; \
- if (jp != 0) { \
- MKL_DC_SWAP(dc_access(a, lda, 0, 0), dc_access(a, lda, 0, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 1, 0), dc_access(a, lda, 1, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 2, 0), dc_access(a, lda, 2, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 3, 0), dc_access(a, lda, 3, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 4, 0), dc_access(a, lda, 4, jp)); \
- MKL_DC_SWAP(dc_access(a, lda, 5, 0), dc_access(a, lda, 5, jp)); \
- } \
- } while (0)
- #define mkl_dc_getri_inverse_upper(n, a, lda, dc_access) \
- do { \
- MKL_INT i, j, ii, jj; \
- mkl_dc_type ajj; \
- for (j = 0; j < n; ++j) { \
- MKL_DC_INV(dc_access(a, lda, j, j), dc_access(a, lda, j, j)); \
- MKL_DC_NEG(ajj, dc_access(a, lda, j, j)); \
- for (jj = 0; jj < j; ++jj) { \
- mkl_dc_type tmp = dc_access(a, lda, jj, j); \
- for (ii = 0; ii < jj; ++ii) { \
- MKL_DC_MUL_ADD(dc_access(a, lda, ii, j), tmp, dc_access(a, lda, ii, jj), dc_access(a, lda, ii, j)); \
- } \
- MKL_DC_MUL(dc_access(a, lda, jj, j), dc_access(a, lda, jj, j), dc_access(a, lda, jj, jj)); \
- } \
- for (i = 0; i < j; ++i) { \
- MKL_DC_MUL(dc_access(a, lda, i, j), dc_access(a, lda, i, j), ajj); \
- } \
- } \
- } while (0)
- #define mkl_dc_getri_solve_lower(n, a, lda, work, dc_access) \
- do { \
- MKL_INT i, j, ii, jj; \
- mkl_dc_type tmp; \
- for (j = n - 1; j >= 0; j--) { \
- for (i = j + 1; i < n; ++i) { \
- work[i] = dc_access(a, lda, i, j); \
- MKL_DC_SET_ZERO(dc_access(a, lda, i, j)); \
- } \
- for (jj = j + 1; jj < n; ++jj) { \
- MKL_DC_NEG(tmp, work[jj]); \
- for (ii = 0; ii < n; ++ii) { \
- MKL_DC_MUL_ADD(dc_access(a, lda, ii, j), tmp, dc_access(a, lda, ii, jj), dc_access(a, lda, ii, j)); \
- } \
- } \
- } \
- } while (0)
- #define mkl_dc_getri_n(n, a, lda, ipiv, work, lwork, info, dc_access) \
- do { \
- MKL_INT i, j, jp; \
- *info = 0; \
- for (i = 0; i < n; ++i) { \
- if (MKL_DC_IS_ZERO(dc_access(a, lda, i, i))) { \
- *info = i + 1; \
- break; \
- } \
- } \
- if (*info == 0) { \
- mkl_dc_getri_inverse_upper(n, a, lda, dc_access); \
- mkl_dc_getri_solve_lower(n, a, lda, work, dc_access); \
- for (j = n - 2; j >= 0; j--) { \
- jp = ipiv[j] - 1; \
- if (jp != j) { \
- for (i = 0; i < n; ++i) { \
- MKL_DC_SWAP(dc_access(a, lda, i, j), dc_access(a, lda, i, jp)); \
- } \
- } \
- } \
- } \
- } while (0)
- #define mkl_dc_getri_generic(n, a, lda, ipiv, work, lwork, info, dc_access) \
- do { \
- if (lwork != -1) { \
- switch (n) { \
- case 0: \
- break; \
- case 1: \
- mkl_dc_getri_1x1(1, a, lda, ipiv, info, dc_access); \
- break; \
- case 2: \
- mkl_dc_getri_2x2(2, a, lda, ipiv, info, dc_access); \
- break; \
- case 3: \
- mkl_dc_getri_3x3(3, a, lda, ipiv, info, dc_access); \
- break; \
- case 4: \
- mkl_dc_getri_4x4(4, a, lda, ipiv, info, dc_access); \
- break; \
- case 5: \
- mkl_dc_getri_5x5(5, a, lda, ipiv, info, dc_access); \
- break; \
- case 6: \
- mkl_dc_getri_6x6(6, a, lda, ipiv, info, dc_access); \
- break; \
- default: \
- mkl_dc_getri_n(n, a, lda, ipiv, work, lwork, info, dc_access); \
- break; \
- } \
- } \
- MKL_DC_CONVERT_INT(work[0], MKL_DC_MAX(1, n)); \
- } while (0)
- static __inline void mkl_dc_getri(MKL_INT n, mkl_dc_type* a, MKL_INT lda, const MKL_INT* ipiv,
- mkl_dc_type* work, MKL_INT lwork, MKL_INT* info)
- {
- mkl_dc_getri_generic(n, a, lda, ipiv, work, lwork, info, MKL_DC_MN);
- }
- #ifndef MKL_DIRECT_CALL_LAPACKE_DISABLE
- static __inline lapack_int mkl_dc_lapacke_getri_convert(int matrix_layout, lapack_int n, mkl_dc_type* a,
- lapack_int lda, const lapack_int* ipiv,
- mkl_dc_type* work, lapack_int lwork)
- {
- lapack_int info = 0;
- if (MKL_DC_GETRI_CHECKSIZE(n)) {
- mkl_dc_type work_small[10];
- lapack_int lwork_small = 10;
- if (!work) {
- work = work_small;
- lwork = lwork_small;
- }
- if (matrix_layout == LAPACK_ROW_MAJOR) {
- mkl_dc_getri_generic(n, a, lda, ipiv, work, lwork, &info, MKL_DC_MT);
- } else {
- mkl_dc_getri_generic(n, a, lda, ipiv, work, lwork, &info, MKL_DC_MN);
- }
- return info;
- }
- if (work) {
- return MKL_DC_CONCAT3(LAPACKE_, MKL_DC_PREC_LETTER, getri_work)(matrix_layout, n, a, lda, ipiv, work, lwork);
- }
- return MKL_DC_CONCAT3(LAPACKE_, MKL_DC_PREC_LETTER, getri)(matrix_layout, n, a, lda, ipiv);
- }
- #endif
- /* ?GEQRF */
- #if defined(MKL_SINGLE) || defined(MKL_DOUBLE)
- #define mkl_dc_geqrf_hh(n, alpha, x, ldx, tau, dc_access) \
- do { \
- MKL_INT ii; \
- MKL_INT n1 = n - 1; \
- mkl_dc_real_type xnorm; \
- mkl_dc_real_type ax; \
- if (n <= 1) { \
- MKL_DC_SET_ZERO(tau); \
- break; \
- } \
- MKL_DC_NRM2_VEC(xnorm, n1, x, ldx, dc_access); \
- if (MKL_DC_IS_R_ZERO(xnorm)) { \
- MKL_DC_SET_ZERO(tau); \
- } else { \
- MKL_INT knt = 0; \
- mkl_dc_type beta; \
- mkl_dc_real_type safmin; \
- mkl_dc_real_type ab; \
- MKL_DC_R_PY2(ax, alpha, xnorm); \
- MKL_DC_R_SIGN(beta, ax, alpha); \
- MKL_DC_NEG(beta, beta); \
- if (!MKL_DC_UNSAFE) { \
- safmin = MKL_DC_XLAMCH("s") / MKL_DC_XLAMCH("e"); \
- if (MKL_DC_ABS(beta) < safmin) { \
- mkl_dc_real_type rsafmn; \
- MKL_DC_INV(rsafmn, safmin); \
- do { \
- knt++; \
- for (ii = 0; ii < n - 1; ++ii) { \
- MKL_DC_MUL(dc_access(x, ldx, ii, 0), dc_access(x, ldx, ii, 0), rsafmn); \
- } \
- MKL_DC_MUL(beta, beta, rsafmn); \
- MKL_DC_MUL(alpha, alpha, rsafmn); \
- } while (MKL_DC_ABS(beta) < safmin); \
- MKL_DC_NRM2_VEC(xnorm, n1, x, ldx, dc_access); \
- MKL_DC_R_PY2(ax, alpha, xnorm); \
- MKL_DC_R_SIGN(beta, ax, alpha); \
- MKL_DC_R_NEG(beta, beta); \
- } \
- } \
- MKL_DC_SUB(tau, beta, alpha); \
- MKL_DC_DIV(tau, tau, beta); \
- MKL_DC_SUB(ab, alpha, beta); \
- for (ii = 0; ii < n - 1; ++ii) { \
- MKL_DC_DIV(dc_access(x, ldx, ii, 0), dc_access(x, ldx, ii, 0), ab); \
- } \
- if (!MKL_DC_UNSAFE) { \
- for (ii = 0; ii < knt; ++ii) { \
- MKL_DC_R_MUL(beta, beta, safmin); \
- } \
- } \
- alpha = beta; \
- } \
- } while (0)
- #endif
- #if defined(MKL_COMPLEX) || defined(MKL_COMPLEX16)
- #define mkl_dc_geqrf_hh(n, alpha, x, ldx, tau, dc_access) \
- do { \
- MKL_INT ii; \
- MKL_INT n1 = n - 1; \
- mkl_dc_real_type xnorm; \
- mkl_dc_real_type ax; \
- mkl_dc_real_type alphr; \
- mkl_dc_real_type alphi; \
- mkl_dc_real_type beta; \
- mkl_dc_real_type rzero; \
- mkl_dc_type temp; \
- mkl_dc_real_type re_temp; \
- mkl_dc_real_type im_temp; \
- MKL_DC_SET_R_ZERO(rzero); \
- if (n <= 0) { \
- MKL_DC_SET_ZERO(tau); \
- break; \
- } \
- MKL_DC_NRM2_VEC(xnorm, n1, x, ldx, dc_access); \
- MKL_DC_REIM(alphr, alphi, alpha); \
- if (MKL_DC_IS_R_ZERO(xnorm) && MKL_DC_IS_R_ZERO(alphi)) { \
- MKL_DC_SET_ZERO(tau); \
- } else { \
- MKL_INT knt = 0; \
- mkl_dc_real_type safmin; \
- MKL_DC_R_PY3(ax, alphr, alphi, xnorm); \
- MKL_DC_R_SIGN(beta, ax, alphr); \
- MKL_DC_R_NEG(beta, beta); \
- if (!MKL_DC_UNSAFE) { \
- safmin = MKL_DC_XLAMCH("s") / MKL_DC_XLAMCH("e"); \
- if (MKL_DC_ABS(beta) < safmin) { \
- mkl_dc_real_type rsafmn; \
- mkl_dc_type c_rsafmn; \
- MKL_DC_R_INV(rsafmn, safmin); \
- MKL_DC_CMPLX(c_rsafmn, rsafmn, rzero); \
- do { \
- knt++; \
- for (ii = 0; ii < n - 1; ++ii) { \
- MKL_DC_MUL(dc_access(x, ldx, ii, 0), dc_access(x, ldx, ii, 0), c_rsafmn); \
- } \
- MKL_DC_R_MUL(beta, beta, rsafmn); \
- MKL_DC_R_MUL(alphi, alphi, rsafmn); \
- MKL_DC_R_MUL(alphr, alphr, rsafmn); \
- } while (MKL_DC_ABS(beta) < safmin); \
- MKL_DC_NRM2_VEC(xnorm, n1, x, ldx, dc_access); \
- MKL_DC_CMPLX(alpha, alphr, alphi); \
- MKL_DC_R_PY3(ax, alphr, alphi, xnorm); \
- MKL_DC_R_SIGN(beta, ax, alphr); \
- MKL_DC_R_NEG(beta, beta); \
- } \
- } \
- MKL_DC_R_SUB(re_temp, beta, alphr); \
- MKL_DC_R_DIV(re_temp, re_temp, beta); \
- MKL_DC_R_NEG(im_temp, alphi); \
- MKL_DC_R_DIV(im_temp, im_temp, beta); \
- MKL_DC_CMPLX(tau, re_temp, im_temp); \
- MKL_DC_CMPLX(temp, beta, rzero); \
- MKL_DC_SUB(temp, alpha, temp); \
- MKL_DC_INV(alpha, temp); \
- for (ii = 0; ii < n - 1; ++ii) { \
- MKL_DC_MUL(dc_access(x, ldx, ii, 0), dc_access(x, ldx, ii, 0), alpha); \
- } \
- if (!MKL_DC_UNSAFE) { \
- for (ii = 0; ii < knt; ++ii) { \
- MKL_DC_R_MUL(beta, beta, safmin); \
- } \
- } \
- MKL_DC_CMPLX(alpha, beta, rzero); \
- } \
- } while (0)
- #endif
- #define mkl_dc_geqrf_apply(m, n, v, ldv, tau, c, ldc, work, dc_access) \
- do { \
- MKL_INT ii; \
- MKL_INT jj; \
- mkl_dc_type mtau; \
- MKL_DC_NEG(mtau, tau); \
- for (ii = 0; ii < n; ++ii) { \
- MKL_DC_SET_ZERO(MKL_DC_ACCESS1D(work, ii)); \
- for (jj = 0; jj < m; ++jj) { \
- mkl_dc_type cji; \
- MKL_DC_CONJ(cji, dc_access(c, ldc, jj, ii)); \
- MKL_DC_MUL_ADD(MKL_DC_ACCESS1D(work, ii), cji, dc_access(v, ldv, jj, 0), MKL_DC_ACCESS1D(work, ii)); \
- } \
- } \
- for (jj = 0; jj < n; ++jj) { \
- mkl_dc_type temp; \
- mkl_dc_type workj; \
- MKL_DC_CONJ(workj, MKL_DC_ACCESS1D(work, jj)); \
- MKL_DC_MUL(temp, mtau, workj); \
- for (ii = 0; ii < m; ++ii) { \
- MKL_DC_MUL_ADD(dc_access(c, ldc, ii, jj), dc_access(v, ldv, ii, 0), temp, dc_access(c, ldc, ii, jj)); \
- } \
- } \
- } while (0)
- #define mkl_dc_geqrf_mxn(m, n, a, lda, tau, work, lwork, info, dc_access) \
- do { \
- MKL_INT i; \
- MKL_INT k = MKL_DC_MIN(m, n); \
- if (lwork == -1) { \
- MKL_DC_CONVERT_INT(MKL_DC_ACCESS1D(work, 0), MKL_DC_MAX(1, n)); \
- } \
- for (i = 0; i < k; ++i) { \
- MKL_INT i1 = MKL_DC_MIN(i + 1, m - 1); \
- MKL_INT mi = m - i; \
- MKL_INT ni1 = n - i - 1; \
- mkl_dc_geqrf_hh((m - i), (dc_access(a, lda, i, i)), (&dc_access(a, lda, i1, i)), lda, MKL_DC_ACCESS1D(tau, i), dc_access); \
- if (i + 1 < n) { \
- mkl_dc_type aii = dc_access(a, lda, i, i); \
- mkl_dc_type taui; \
- MKL_DC_SET_ONE(dc_access(a, lda, i, i)); \
- MKL_DC_CONJ(taui, MKL_DC_ACCESS1D(tau, i)); \
- mkl_dc_geqrf_apply(mi, ni1, (&dc_access(a, lda, i, i)), lda, taui, (&dc_access(a, lda, i, i + 1)), lda, work, dc_access); \
- dc_access(a, lda, i, i) = aii; \
- } \
- } \
- MKL_DC_CONVERT_INT(MKL_DC_ACCESS1D(work, 0), MKL_DC_MAX(1, n)); \
- } while (0)
- #define mkl_dc_geqrf_generic(m, n, a, lda, tau, work, lwork, info, dc_access) \
- do { \
- *info = 0; \
- if (lwork != -1) { \
- if (m == n && m <= 5) { \
- switch (n) { \
- case 0: \
- break; \
- case 1: \
- mkl_dc_geqrf_mxn(1, 1, a, lda, tau, work, lwork, info, dc_access); \
- break; \
- case 2: \
- mkl_dc_geqrf_mxn(2, 2, a, lda, tau, work, lwork, info, dc_access); \
- break; \
- case 3: \
- mkl_dc_geqrf_mxn(3, 3, a, lda, tau, work, lwork, info, dc_access); \
- break; \
- case 4: \
- mkl_dc_geqrf_mxn(4, 4, a, lda, tau, work, lwork, info, dc_access); \
- break; \
- case 5: \
- mkl_dc_geqrf_mxn(5, 5, a, lda, tau, work, lwork, info, dc_access); \
- break; \
- } \
- } else { \
- mkl_dc_geqrf_mxn(m, n, a, lda, tau, work, lwork, info, dc_access); \
- } \
- } else { \
- MKL_DC_CONVERT_INT(MKL_DC_ACCESS1D(work, 0), MKL_DC_MAX(1, n)); \
- } \
- } while (0)
- static __inline void mkl_dc_geqrf(MKL_INT m, MKL_INT n, mkl_dc_type* a, MKL_INT lda, mkl_dc_type* tau,
- mkl_dc_type* work, MKL_INT lwork, MKL_INT* info)
- {
- mkl_dc_geqrf_generic(m, n, a, lda, tau, work, lwork, info, MKL_DC_MN);
- }
- #ifndef MKL_DIRECT_CALL_LAPACKE_DISABLE
- static __inline lapack_int mkl_dc_lapacke_geqrf_convert(int matrix_layout, lapack_int m, lapack_int n,
- mkl_dc_type* a, lapack_int lda,
- mkl_dc_type* tau, mkl_dc_type* work, lapack_int lwork)
- {
- lapack_int info = 0;
- if (MKL_DC_GEQRF_CHECKSIZE(m, n)) {
- mkl_dc_type work_small[10];
- lapack_int lwork_small = 10;
- if (!work) {
- work = work_small;
- lwork = lwork_small;
- }
- if (matrix_layout == LAPACK_ROW_MAJOR) {
- mkl_dc_geqrf_generic(m, n, a, lda, tau, work, lwork, &info, MKL_DC_MT);
- } else {
- mkl_dc_geqrf_generic(m, n, a, lda, tau, work, lwork, &info, MKL_DC_MN);
- }
- return info;
- }
- if (work) {
- return MKL_DC_CONCAT3(LAPACKE_, MKL_DC_PREC_LETTER, geqrf_work)(matrix_layout, m, n, a, lda, tau, work, lwork);
- }
- return MKL_DC_CONCAT3(LAPACKE_, MKL_DC_PREC_LETTER, geqrf)(matrix_layout, m, n, a, lda, tau);
- }
- #endif
- /* ?POTRF */
- #define mkl_dc_potrf_n(uplo, n, a, lda, info, dc_access) \
- do { \
- MKL_INT i, j, ii, jj; \
- if (MKL_DC_MisU(uplo)) { \
- for (j = 0; j < n; ++j) { \
- mkl_dc_type ajj = dc_access(a, lda, j, j); \
- MKL_DC_ZERO_IMAG(ajj); \
- for (i = 0; i < j; ++i) { \
- MKL_DC_SUB_MUL_CONJ(ajj, ajj, dc_access(a, lda, i, j), dc_access(a, lda, i, j)); \
- } \
- if (MKL_DC_NON_POS(ajj)) { \
- dc_access(a, lda, j, j) = ajj; \
- *info = j + 1; \
- break; \
- } \
- MKL_DC_SQRT(ajj, ajj); \
- dc_access(a, lda, j, j) = ajj; \
- MKL_DC_INV(ajj, ajj); \
- if (j < n - 1) { \
- for (jj = j + 1; jj < n; ++jj) { \
- mkl_dc_type aj_jj = dc_access(a, lda, j, jj); \
- for (ii = 0; ii < j; ++ii) { \
- MKL_DC_SUB_MUL_CONJ(aj_jj, aj_jj, \
- dc_access(a, lda, ii, jj), dc_access(a, lda, ii, j)); \
- } \
- dc_access(a, lda, j, jj) = aj_jj; \
- } \
- for (i = j + 1; i < n; ++i) { \
- MKL_DC_MUL(dc_access(a, lda, j, i), dc_access(a, lda, j, i), ajj); \
- } \
- } \
- } \
- } else { \
- for (j = 0; j < n; ++j) { \
- mkl_dc_type ajj = dc_access(a, lda, j, j); \
- MKL_DC_ZERO_IMAG(ajj); \
- for (i = 0; i < j; ++i) { \
- MKL_DC_SUB_MUL_CONJ(ajj, ajj, dc_access(a, lda, j, i), dc_access(a, lda, j, i)); \
- } \
- if (MKL_DC_NON_POS(ajj)) { \
- dc_access(a, lda, j, j) = ajj; \
- *info = j + 1; \
- break; \
- } \
- MKL_DC_SQRT(ajj, ajj); \
- dc_access(a, lda, j, j) = ajj; \
- MKL_DC_INV(ajj, ajj); \
- if (j < n - 1) { \
- for (ii = 0; ii < j; ++ii) { \
- for (jj = j + 1; jj < n; ++jj) { \
- MKL_DC_SUB_MUL_CONJ(dc_access(a, lda, jj, j), dc_access(a, lda, jj, j), \
- dc_access(a, lda, jj, ii), dc_access(a, lda, j, ii)); \
- } \
- } \
- for (i = j + 1; i < n; ++i) { \
- MKL_DC_MUL(dc_access(a, lda, i, j), dc_access(a, lda, i, j), ajj); \
- } \
- } \
- } \
- } \
- } while (0) \
- #define mkl_dc_potrf_1(uplo, n, a, lda, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- MKL_DC_ZERO_IMAG(a00); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- } while (0) \
- #define mkl_dc_potrf_2(uplo, n, a, lda, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- MKL_DC_ZERO_IMAG(a00); \
- MKL_DC_ZERO_IMAG(a11); \
- if (MKL_DC_MisU(uplo)) { \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a01, a01, a00); \
- dc_access(a, lda, 0, 1) = a01; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a01, a01); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- } else { \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a10, a10, a00); \
- dc_access(a, lda, 1, 0) = a10; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a10, a10); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- } \
- } while (0)
- #define mkl_dc_potrf_3(uplo, n, a, lda, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- MKL_DC_ZERO_IMAG(a00); \
- MKL_DC_ZERO_IMAG(a11); \
- MKL_DC_ZERO_IMAG(a22); \
- if (MKL_DC_MisU(uplo)) { \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a01, a01, a00); \
- dc_access(a, lda, 0, 1) = a01; \
- MKL_DC_MUL(a02, a02, a00); \
- dc_access(a, lda, 0, 2) = a02; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a01, a01); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- MKL_DC_SUB_MUL_CONJ(a12, a12, a02, a01); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_MUL(a12, a12, a11); \
- dc_access(a, lda, 1, 2) = a12; \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a02, a02); \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a12, a12); \
- if (MKL_DC_NON_POS(a22)) { \
- dc_access(a, lda, 2, 2) = a22; \
- *info = 3; \
- break; \
- } \
- MKL_DC_SQRT(a22, a22); \
- dc_access(a, lda, 2, 2) = a22; \
- } else { \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a10, a10, a00); \
- dc_access(a, lda, 1, 0) = a10; \
- MKL_DC_MUL(a20, a20, a00); \
- dc_access(a, lda, 2, 0) = a20; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a10, a10); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- MKL_DC_SUB_MUL_CONJ(a21, a21, a20, a10); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_MUL(a21, a21, a11); \
- dc_access(a, lda, 2, 1) = a21; \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a20, a20); \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a21, a21); \
- if (MKL_DC_NON_POS(a22)) { \
- dc_access(a, lda, 2, 2) = a22; \
- *info = 3; \
- break; \
- } \
- MKL_DC_SQRT(a22, a22); \
- dc_access(a, lda, 2, 2) = a22; \
- } \
- } while (0)
- #define mkl_dc_potrf_4(uplo, n, a, lda, info, dc_access) \
- do { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- MKL_DC_ZERO_IMAG(a00); \
- MKL_DC_ZERO_IMAG(a11); \
- MKL_DC_ZERO_IMAG(a22); \
- MKL_DC_ZERO_IMAG(a33); \
- if (MKL_DC_MisU(uplo)) { \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a01, a01, a00); \
- dc_access(a, lda, 0, 1) = a01; \
- MKL_DC_MUL(a02, a02, a00); \
- dc_access(a, lda, 0, 2) = a02; \
- MKL_DC_MUL(a03, a03, a00); \
- dc_access(a, lda, 0, 3) = a03; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a01, a01); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- MKL_DC_SUB_MUL_CONJ(a12, a12, a02, a01); \
- MKL_DC_SUB_MUL_CONJ(a13, a13, a03, a01); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_MUL(a12, a12, a11); \
- dc_access(a, lda, 1, 2) = a12; \
- MKL_DC_MUL(a13, a13, a11); \
- dc_access(a, lda, 1, 3) = a13; \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a02, a02); \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a12, a12); \
- if (MKL_DC_NON_POS(a22)) { \
- dc_access(a, lda, 2, 2) = a22; \
- *info = 3; \
- break; \
- } \
- MKL_DC_SQRT(a22, a22); \
- dc_access(a, lda, 2, 2) = a22; \
- MKL_DC_SUB_MUL_CONJ(a23, a23, a03, a02); \
- MKL_DC_SUB_MUL_CONJ(a23, a23, a13, a12); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_MUL(a23, a23, a22); \
- dc_access(a, lda, 2, 3) = a23; \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a03, a03); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a13, a13); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a23, a23); \
- if (MKL_DC_NON_POS(a33)) { \
- dc_access(a, lda, 3, 3) = a33; \
- *info = 4; \
- break; \
- } \
- MKL_DC_SQRT(a33, a33); \
- dc_access(a, lda, 3, 3) = a33; \
- } else { \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a10, a10, a00); \
- dc_access(a, lda, 1, 0) = a10; \
- MKL_DC_MUL(a20, a20, a00); \
- dc_access(a, lda, 2, 0) = a20; \
- MKL_DC_MUL(a30, a30, a00); \
- dc_access(a, lda, 3, 0) = a30; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a10, a10); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- MKL_DC_SUB_MUL_CONJ(a21, a21, a20, a10); \
- MKL_DC_SUB_MUL_CONJ(a31, a31, a30, a10); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_MUL(a21, a21, a11); \
- dc_access(a, lda, 2, 1) = a21; \
- MKL_DC_MUL(a31, a31, a11); \
- dc_access(a, lda, 3, 1) = a31; \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a20, a20); \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a21, a21); \
- if (MKL_DC_NON_POS(a22)) { \
- dc_access(a, lda, 2, 2) = a22; \
- *info = 3; \
- break; \
- } \
- MKL_DC_SQRT(a22, a22); \
- dc_access(a, lda, 2, 2) = a22; \
- MKL_DC_SUB_MUL_CONJ(a32, a32, a30, a20); \
- MKL_DC_SUB_MUL_CONJ(a32, a32, a31, a21); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_MUL(a32, a32, a22); \
- dc_access(a, lda, 3, 2) = a32; \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a30, a30); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a31, a31); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a32, a32); \
- if (MKL_DC_NON_POS(a33)) { \
- dc_access(a, lda, 3, 3) = a33; \
- *info = 4; \
- break; \
- } \
- MKL_DC_SQRT(a33, a33); \
- dc_access(a, lda, 3, 3) = a33; \
- } \
- } while (0)
- #define mkl_dc_potrf_5(uplo, n, a, lda, info, dc_access) \
- { \
- mkl_dc_type a00 = dc_access(a, lda, 0, 0); \
- mkl_dc_type a11 = dc_access(a, lda, 1, 1); \
- mkl_dc_type a22 = dc_access(a, lda, 2, 2); \
- mkl_dc_type a33 = dc_access(a, lda, 3, 3); \
- mkl_dc_type a44 = dc_access(a, lda, 4, 4); \
- MKL_DC_ZERO_IMAG(a00); \
- MKL_DC_ZERO_IMAG(a11); \
- MKL_DC_ZERO_IMAG(a22); \
- MKL_DC_ZERO_IMAG(a33); \
- MKL_DC_ZERO_IMAG(a44); \
- if (MKL_DC_MisU(uplo)) { \
- mkl_dc_type a01 = dc_access(a, lda, 0, 1); \
- mkl_dc_type a02 = dc_access(a, lda, 0, 2); \
- mkl_dc_type a12 = dc_access(a, lda, 1, 2); \
- mkl_dc_type a03 = dc_access(a, lda, 0, 3); \
- mkl_dc_type a13 = dc_access(a, lda, 1, 3); \
- mkl_dc_type a23 = dc_access(a, lda, 2, 3); \
- mkl_dc_type a04 = dc_access(a, lda, 0, 4); \
- mkl_dc_type a14 = dc_access(a, lda, 1, 4); \
- mkl_dc_type a24 = dc_access(a, lda, 2, 4); \
- mkl_dc_type a34 = dc_access(a, lda, 3, 4); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a01, a01, a00); \
- MKL_DC_MUL(a02, a02, a00); \
- MKL_DC_MUL(a03, a03, a00); \
- MKL_DC_MUL(a04, a04, a00); \
- dc_access(a, lda, 0, 1) = a01; \
- dc_access(a, lda, 0, 2) = a02; \
- dc_access(a, lda, 0, 3) = a03; \
- dc_access(a, lda, 0, 4) = a04; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a01, a01); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- MKL_DC_SUB_MUL_CONJ(a12, a12, a02, a01); \
- MKL_DC_SUB_MUL_CONJ(a13, a13, a03, a01); \
- MKL_DC_SUB_MUL_CONJ(a14, a14, a04, a01); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_MUL(a12, a12, a11); \
- MKL_DC_MUL(a13, a13, a11); \
- MKL_DC_MUL(a14, a14, a11); \
- dc_access(a, lda, 1, 2) = a12; \
- dc_access(a, lda, 1, 3) = a13; \
- dc_access(a, lda, 1, 4) = a14; \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a02, a02); \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a12, a12); \
- if (MKL_DC_NON_POS(a22)) { \
- dc_access(a, lda, 2, 2) = a22; \
- *info = 3; \
- break; \
- } \
- MKL_DC_SQRT(a22, a22); \
- dc_access(a, lda, 2, 2) = a22; \
- MKL_DC_SUB_MUL_CONJ(a23, a23, a03, a02); \
- MKL_DC_SUB_MUL_CONJ(a23, a23, a13, a12); \
- MKL_DC_SUB_MUL_CONJ(a24, a24, a04, a02); \
- MKL_DC_SUB_MUL_CONJ(a24, a24, a14, a12); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_MUL(a23, a23, a22); \
- MKL_DC_MUL(a24, a24, a22); \
- dc_access(a, lda, 2, 3) = a23; \
- dc_access(a, lda, 2, 4) = a24; \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a03, a03); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a13, a13); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a23, a23); \
- if (MKL_DC_NON_POS(a33)) { \
- dc_access(a, lda, 3, 3) = a33; \
- *info = 4; \
- break; \
- } \
- MKL_DC_SQRT(a33, a33); \
- dc_access(a, lda, 3, 3) = a33; \
- MKL_DC_SUB_MUL_CONJ(a34, a34, a04, a03); \
- MKL_DC_SUB_MUL_CONJ(a34, a34, a14, a13); \
- MKL_DC_SUB_MUL_CONJ(a34, a34, a24, a23); \
- MKL_DC_INV(a33, a33); \
- MKL_DC_MUL(a34, a34, a33); \
- dc_access(a, lda, 3, 4) = a34; \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a04, a04); \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a14, a14); \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a24, a24); \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a34, a34); \
- if (MKL_DC_NON_POS(a44)) { \
- dc_access(a, lda, 4, 4) = a44; \
- *info = 5; \
- break; \
- } \
- MKL_DC_SQRT(a44, a44); \
- dc_access(a, lda, 4, 4) = a44; \
- } else { \
- mkl_dc_type a10 = dc_access(a, lda, 1, 0); \
- mkl_dc_type a20 = dc_access(a, lda, 2, 0); \
- mkl_dc_type a30 = dc_access(a, lda, 3, 0); \
- mkl_dc_type a40 = dc_access(a, lda, 4, 0); \
- mkl_dc_type a21 = dc_access(a, lda, 2, 1); \
- mkl_dc_type a43 = dc_access(a, lda, 4, 3); \
- mkl_dc_type a31 = dc_access(a, lda, 3, 1); \
- mkl_dc_type a41 = dc_access(a, lda, 4, 1); \
- mkl_dc_type a32 = dc_access(a, lda, 3, 2); \
- mkl_dc_type a42 = dc_access(a, lda, 4, 2); \
- if (MKL_DC_NON_POS(a00)) { \
- dc_access(a, lda, 0, 0) = a00; \
- *info = 1; \
- break; \
- } \
- MKL_DC_SQRT(a00, a00); \
- dc_access(a, lda, 0, 0) = a00; \
- MKL_DC_INV(a00, a00); \
- MKL_DC_MUL(a10, a10, a00); \
- MKL_DC_MUL(a20, a20, a00); \
- MKL_DC_MUL(a30, a30, a00); \
- MKL_DC_MUL(a40, a40, a00); \
- dc_access(a, lda, 1, 0) = a10; \
- dc_access(a, lda, 2, 0) = a20; \
- dc_access(a, lda, 3, 0) = a30; \
- dc_access(a, lda, 4, 0) = a40; \
- MKL_DC_SUB_MUL_CONJ(a11, a11, a10, a10); \
- if (MKL_DC_NON_POS(a11)) { \
- dc_access(a, lda, 1, 1) = a11; \
- *info = 2; \
- break; \
- } \
- MKL_DC_SQRT(a11, a11); \
- dc_access(a, lda, 1, 1) = a11; \
- MKL_DC_SUB_MUL_CONJ(a21, a21, a20, a10); \
- MKL_DC_SUB_MUL_CONJ(a31, a31, a30, a10); \
- MKL_DC_SUB_MUL_CONJ(a41, a41, a40, a10); \
- MKL_DC_INV(a11, a11); \
- MKL_DC_MUL(a21, a21, a11); \
- MKL_DC_MUL(a31, a31, a11); \
- MKL_DC_MUL(a41, a41, a11); \
- dc_access(a, lda, 2, 1) = a21; \
- dc_access(a, lda, 3, 1) = a31; \
- dc_access(a, lda, 4, 1) = a41; \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a20, a20); \
- MKL_DC_SUB_MUL_CONJ(a22, a22, a21, a21); \
- if (MKL_DC_NON_POS(a22)) { \
- dc_access(a, lda, 2, 2) = a22; \
- *info = 3; \
- break; \
- } \
- MKL_DC_SQRT(a22, a22); \
- dc_access(a, lda, 2, 2) = a22; \
- MKL_DC_SUB_MUL_CONJ(a32, a32, a30, a20); \
- MKL_DC_SUB_MUL_CONJ(a32, a32, a31, a21); \
- MKL_DC_SUB_MUL_CONJ(a42, a42, a40, a20); \
- MKL_DC_SUB_MUL_CONJ(a42, a42, a41, a21); \
- MKL_DC_INV(a22, a22); \
- MKL_DC_MUL(a32, a32, a22); \
- MKL_DC_MUL(a42, a42, a22); \
- dc_access(a, lda, 3, 2) = a32; \
- dc_access(a, lda, 4, 2) = a42; \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a30, a30); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a31, a31); \
- MKL_DC_SUB_MUL_CONJ(a33, a33, a32, a32); \
- if (MKL_DC_NON_POS(a33)) { \
- dc_access(a, lda, 3, 3) = a33; \
- *info = 4; \
- break; \
- } \
- MKL_DC_SQRT(a33, a33); \
- dc_access(a, lda, 3, 3) = a33; \
- MKL_DC_SUB_MUL_CONJ(a43, a43, a40, a30); \
- MKL_DC_SUB_MUL_CONJ(a43, a43, a41, a31); \
- MKL_DC_SUB_MUL_CONJ(a43, a43, a42, a32); \
- MKL_DC_INV(a33, a33); \
- MKL_DC_MUL(a43, a43, a33); \
- dc_access(a, lda, 4, 3) = a43; \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a40, a40); \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a41, a41); \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a42, a42); \
- MKL_DC_SUB_MUL_CONJ(a44, a44, a43, a43); \
- if (MKL_DC_NON_POS(a44)) { \
- dc_access(a, lda, 4, 4) = a44; \
- *info = 5; \
- break; \
- } \
- MKL_DC_SQRT(a44, a44); \
- dc_access(a, lda, 4, 4) = a44; \
- } \
- } while (0)
- #define mkl_dc_potrf_uplo(uplo, n, a, lda, info, dc_access) \
- do { \
- if (n <= 5) { \
- switch (n) { \
- case 0: \
- break; \
- case 1: \
- mkl_dc_potrf_1(uplo, n, a, lda, info, dc_access); \
- break; \
- case 2: \
- mkl_dc_potrf_2(uplo, n, a, lda, info, dc_access); \
- break; \
- case 3: \
- mkl_dc_potrf_3(uplo, n, a, lda, info, dc_access); \
- break; \
- case 4: \
- mkl_dc_potrf_4(uplo, n, a, lda, info, dc_access); \
- break; \
- case 5: \
- mkl_dc_potrf_5(uplo, n, a, lda, info, dc_access); \
- break; \
- } \
- } else { \
- mkl_dc_potrf_n(uplo, n, a, lda, info, dc_access); \
- } \
- } while (0)
- #define mkl_dc_potrf_generic(uplo, n, a, lda, info, dc_access) \
- do { \
- int uuplo = MKL_DC_MisU(uplo); \
- *info = 0; \
- if (n != 0) { \
- if (uuplo) { \
- mkl_dc_potrf_uplo('U', n, a, lda, info, dc_access); \
- } else { \
- mkl_dc_potrf_uplo('L', n, a, lda, info, dc_access); \
- } \
- } \
- } while (0)
- static __inline void mkl_dc_potrf(char uplo, MKL_INT n, mkl_dc_type* a, MKL_INT lda, MKL_INT* info)
- {
- mkl_dc_potrf_generic(uplo, n, a, lda, info, MKL_DC_MN);
- }
- #ifndef MKL_DIRECT_CALL_LAPACKE_DISABLE
- static __inline lapack_int mkl_dc_lapacke_potrf_convert(int matrix_layout, char uplo, lapack_int n, mkl_dc_type* a, lapack_int lda)
- {
- lapack_int info = 0;
- if (MKL_DC_POTRF_CHECKSIZE(n)) {
- if (matrix_layout == LAPACK_ROW_MAJOR) {
- mkl_dc_potrf_generic(uplo, n, a, lda, &info, MKL_DC_MT);
- } else {
- mkl_dc_potrf_generic(uplo, n, a, lda, &info, MKL_DC_MN);
- }
- return info;
- }
- return MKL_DC_CONCAT3(LAPACKE_, MKL_DC_PREC_LETTER, potrf)(matrix_layout, uplo, n, a, lda);
- }
- #endif
|