Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.
7034 Discussions

?getrf Failing to Detect Singular Matrices

Kevin_T_
Beginner
‎05-14-2024 06:41 PM
3,142 Views

I'm encountering a widespread issue where MKL's getrf routines can fail to detect singular matrices.

Across several machines (different OSs, CPUs, MKL versions, and gfortran versions) INFO returns zero for obviously singular matrices. I've attached a minimal (ish) working example below. I'm compiling (with gfortran) as:

 

 

 

 

gfortran getrfSingularMWE.f90  -m64  -L${MKLROOT}/lib/intel64 -Wl,--no-as-needed -lmkl_gf_lp64 -lmkl_gnu_thread -lmkl_core -lgomp -lpthread -lm -ldl  -o getrfSingularMWE.x

 

 

 

 


Some machines miss more singular matrices than others. I have found a machine where MKL works as intended, but I'm struggling to determine what might be special about it - identical compiler and MKL versions on other machines will fail. The issue is also not confined to GFortran, Intel Fortran compilers also have failures.

Is there anything I might try to diagnose the circumstances under which MKL's getrf will fail?

 

 

program mkl_check_singular

  use iso_fortran_env, only: REAL32, REAL64

  implicit none

!=====================================================================
! This program checks the reporting of singular LU'd matrices using all
! combinations of single/double, real/complex, and array/matrix for an
! 8x12 matrix.
!
! Singular matrices are created by assigning multiple columns to be the
! same (matrix(:, 2) = matrix(:, 1), etc.).
!
! Different LU algorithms can fail to find slightly singular matrices 
! (i.e. only a few identical columns) in different ways, failing to 
! recognize some while recognizing others that alternative algorithms 
! don't find. If an algorithm fails to recognize a singular U when all 
! columns are identical, however, that algorithm should be treated with 
! suspicion. 
!=====================================================================

  integer :: numRows = 12, numColumns = 8, LDA = 12, INFO, &
             linIndx, rowIndx, colIndx, LMajorVer, LMinorVer, LPatchVer, &
             warningStatus, errorStatus 

  integer, dimension(8) :: IPIV
        
  real(kind=REAL32),       dimension(96) :: sArray, sArraySingular
  real(kind=REAL64),       dimension(96) :: dArray, dArraySingular
  complex(kind=REAL32),    dimension(96) :: cArray, cArraySingular
  complex(kind=REAL64),    dimension(96) :: zArray, zArraySingular

  real(kind=REAL32),    dimension(12,  :: sMatrix, sMatrixSingular
  real(kind=REAL64),    dimension(12,  :: dMatrix, dMatrixSingular
  complex(kind=REAL32), dimension(12,  :: cMatrix, cMatrixSingular
  complex(kind=REAL64), dimension(12,  :: zMatrix, zMatrixSingular

!=====================================================================
! Setup Arrays
!=====================================================================

  ! Single Real
  sArray = [ &
    8.984861_REAL32, 1.181552_REAL32, 9.884179_REAL32, 5.399821_REAL32, 7.069174_REAL32, 9.994916_REAL32, & ! 6
    2.878493_REAL32, 4.145225_REAL32, 4.648399_REAL32, 7.639571_REAL32, 8.182040_REAL32, 1.002215_REAL32, & ! 12
    1.781170_REAL32, 3.596349_REAL32, 0.567047_REAL32, 5.218857_REAL32, 3.358490_REAL32, 1.756690_REAL32, & ! 18
    2.089467_REAL32, 9.051536_REAL32, 6.753912_REAL32, 4.684682_REAL32, 9.121325_REAL32, 1.040116_REAL32, & ! 24
    7.455461_REAL32, 7.362675_REAL32, 5.618614_REAL32, 1.841941_REAL32, 5.972114_REAL32, 2.999370_REAL32, & ! 30
    1.341229_REAL32, 2.126015_REAL32, 8.949417_REAL32, 0.714528_REAL32, 2.424866_REAL32, 0.537544_REAL32, & ! 36
    4.417221_REAL32, 0.132832_REAL32, 8.971914_REAL32, 1.966582_REAL32, 0.933705_REAL32, 3.073669_REAL32, & ! 42
    4.560577_REAL32, 1.016694_REAL32, 9.953897_REAL32, 3.320928_REAL32, 2.973468_REAL32, 0.620452_REAL32, & ! 48
    2.982440_REAL32, 0.463513_REAL32, 5.054281_REAL32, 7.614259_REAL32, 6.310700_REAL32, 0.898917_REAL32, & ! 54
    0.808624_REAL32, 7.772405_REAL32, 9.051347_REAL32, 5.337720_REAL32, 1.091542_REAL32, 8.258089_REAL32, & ! 60
    3.380977_REAL32, 2.939731_REAL32, 7.463134_REAL32, 0.103366_REAL32, 0.484473_REAL32, 6.679161_REAL32, & ! 66
    6.034680_REAL32, 5.261025_REAL32, 7.297094_REAL32, 7.072535_REAL32, 7.813771_REAL32, 2.879770_REAL32, & ! 72
    6.925320_REAL32, 5.566698_REAL32, 3.965208_REAL32, 0.615907_REAL32, 7.801755_REAL32, 3.375839_REAL32, & ! 78
    6.078659_REAL32, 7.412540_REAL32, 1.048132_REAL32, 1.278884_REAL32, 5.495401_REAL32, 4.852294_REAL32, & ! 84
    8.904757_REAL32, 7.989603_REAL32, 7.343411_REAL32, 0.513319_REAL32, 0.728853_REAL32, 0.885275_REAL32, & ! 90
    7.983509_REAL32, 9.430081_REAL32, 6.837156_REAL32, 1.320830_REAL32, 7.227245_REAL32, 1.103535_REAL32]   ! 96

  ! Single Complex
  cArray = [ &
    (1.818470_REAL32, 2.638029_REAL32), (1.455390_REAL32, 1.360686_REAL32), (8.692922_REAL32, 5.797046_REAL32), & ! 3
    (5.498602_REAL32, 1.449548_REAL32), (8.530311_REAL32, 6.220551_REAL32), (3.509524_REAL32, 5.132495_REAL32), & ! 6 
    (4.018080_REAL32, 0.759667_REAL32), (2.399162_REAL32, 1.233189_REAL32), (1.839078_REAL32, 2.399525_REAL32), & ! 9
    (4.172671_REAL32, 0.496544_REAL32), (9.027161_REAL32, 9.447872_REAL32), (4.908641_REAL32, 4.892526_REAL32), & ! 12
    (3.377194_REAL32, 9.000538_REAL32), (3.692468_REAL32, 1.112028_REAL32), (7.802521_REAL32, 3.897388_REAL32), & ! 15
    (2.416913_REAL32, 4.039121_REAL32), (0.964545_REAL32, 1.319733_REAL32), (9.420506_REAL32, 9.561345_REAL32), & ! 18
    (5.752086_REAL32, 0.597795_REAL32), (2.347799_REAL32, 3.531586_REAL32), (8.211940_REAL32, 0.154034_REAL32), & ! 21
    (0.430238_REAL32, 1.689900_REAL32), (6.491155_REAL32, 7.317224_REAL32), (6.477460_REAL32, 4.509237_REAL32), & ! 24
    (5.470089_REAL32, 2.963208_REAL32), (7.446928_REAL32, 1.889550_REAL32), (6.867754_REAL32, 1.835112_REAL32), & ! 27
    (3.684846_REAL32, 6.256186_REAL32), (7.802274_REAL32, 0.811258_REAL32), (9.293860_REAL32, 7.757127_REAL32), & ! 30
    (4.867916_REAL32, 4.358586_REAL32), (4.467837_REAL32, 3.063495_REAL32), (5.085087_REAL32, 5.107716_REAL32), & ! 33
    (8.176277_REAL32, 7.948314_REAL32), (6.443181_REAL32, 3.786094_REAL32), (8.115805_REAL32, 5.328256_REAL32), & ! 36
    (3.507271_REAL32, 9.390016_REAL32), (8.759428_REAL32, 5.501563_REAL32), (6.224751_REAL32, 5.870447_REAL32), & ! 39
    (2.077423_REAL32, 3.012463_REAL32), (4.709233_REAL32, 2.304882_REAL32), (8.443088_REAL32, 1.947643_REAL32), & ! 42
    (2.259218_REAL32, 1.707080_REAL32), (2.276643_REAL32, 4.356987_REAL32), (3.111023_REAL32, 9.233796_REAL32), & ! 25
    (4.302074_REAL32, 1.848163_REAL32), (9.048810_REAL32, 9.797484_REAL32), (4.388700_REAL32, 1.111192_REAL32), & ! 48
    (2.580647_REAL32, 4.087198_REAL32), (5.948961_REAL32, 2.622117_REAL32), (6.028431_REAL32, 7.112158_REAL32), & ! 51
    (2.217467_REAL32, 1.174177_REAL32), (2.966759_REAL32, 3.187783_REAL32), (4.241668_REAL32, 5.078583_REAL32), & ! 54
    (0.855158_REAL32, 2.624822_REAL32), (8.010146_REAL32, 0.292203_REAL32), (9.288541_REAL32, 7.303309_REAL32), & ! 57
    (4.886090_REAL32, 5.785251_REAL32), (2.372836_REAL32, 4.588488_REAL32), (9.630885_REAL32, 5.468057_REAL32), & ! 60
    (5.211358_REAL32, 2.315944_REAL32), (4.888977_REAL32, 6.240601_REAL32), (6.791355_REAL32, 3.955152_REAL32), & ! 63
    (3.674366_REAL32, 9.879820_REAL32), (0.377389_REAL32, 8.851680_REAL32), (9.132868_REAL32, 7.961839_REAL32), & ! 66
    (0.987123_REAL32, 2.618712_REAL32), (3.353568_REAL32, 6.797280_REAL32), (1.365531_REAL32, 7.212275_REAL32), & ! 69
    (1.067619_REAL32, 6.537573_REAL32), (4.941739_REAL32, 7.790517_REAL32), (7.150371_REAL32, 9.037206_REAL32), & ! 72
    (8.909225_REAL32, 3.341631_REAL32), (6.987458_REAL32, 1.978098_REAL32), (0.305409_REAL32, 7.440743_REAL32), & ! 75
    (5.000224_REAL32, 4.799221_REAL32), (9.047222_REAL32, 6.098666_REAL32), (6.176664_REAL32, 8.594423_REAL32), & ! 78
    (8.054894_REAL32, 5.767215_REAL32), (1.829225_REAL32, 2.399320_REAL32), (8.865119_REAL32, 0.286742_REAL32), & ! 81
    (4.899014_REAL32, 1.679271_REAL32), (9.786806_REAL32, 7.126945_REAL32), (5.004716_REAL32, 4.710884_REAL32), & ! 84
    (0.596189_REAL32, 6.819719_REAL32), (0.424311_REAL32, 0.714455_REAL32), (5.216498_REAL32, 0.967300_REAL32), & ! 87
    (8.181486_REAL32, 8.175471_REAL32), (7.224396_REAL32, 1.498654_REAL32), (6.596053_REAL32, 5.185949_REAL32), & ! 90
    (9.729746_REAL32, 6.489915_REAL32), (8.003306_REAL32, 4.537977_REAL32), (4.323915_REAL32, 8.253138_REAL32), & ! 93
    (0.834698_REAL32, 1.331710_REAL32), (1.733886_REAL32, 3.909378_REAL32), (8.313797_REAL32, 8.033644_REAL32)]   ! 96

  ! Double Real
  dArray = [ &
    9.053306072443640_REAL64, 6.401942483322220_REAL64, 1.629395714406698_REAL64, 5.659105910016573_REAL64, & ! 4 
    9.316161369771947_REAL64, 7.831018836846311_REAL64, 6.856871369789888_REAL64, 4.662192849827422_REAL64, & ! 8
    2.603180956665296_REAL64, 5.692681710770191_REAL64, 2.487706628880768_REAL64, 3.193015915208970_REAL64, & ! 12 
    9.108022322169406_REAL64, 8.852200409944105_REAL64, 7.945894241658521_REAL64, 9.258098071668771_REAL64, & ! 16 
    1.788399507758780_REAL64, 5.175411482364652_REAL64, 6.270053164093095_REAL64, 9.131823459009045_REAL64, & ! 20 
    6.639682538388079_REAL64, 3.891928229926131_REAL64, 7.400075814161521_REAL64, 8.176348597087308_REAL64, & ! 24 
    6.003447931764931_REAL64, 0.849970645027028_REAL64, 9.223580010900541_REAL64, 0.535978153165582_REAL64, & ! 28
    5.270249526773659_REAL64, 1.188532766197476_REAL64, 3.801429902410049_REAL64, 8.128325330151092_REAL64, & ! 32 
    2.440959155716900_REAL64, 8.844226590935749_REAL64, 7.126467947591497_REAL64, 3.781484380821433_REAL64, & ! 36 
    2.489196109959332_REAL64, 2.528537452392730_REAL64, 7.672435945269350_REAL64, 0.498618850472403_REAL64, & ! 40 
    6.852885676089091_REAL64, 6.202781221510013_REAL64, 7.466846299687251_REAL64, 9.772556483456693_REAL64, & ! 44 
    3.839135126606154_REAL64, 2.602056023137224_REAL64, 8.774695360141239_REAL64, 8.060959835106466_REAL64, & ! 48 
    4.611210797796994_REAL64, 0.909616879242212_REAL64, 5.642688252731981_REAL64, 1.873828891486068_REAL64, & ! 52 
    5.316895858947326_REAL64, 3.550333217733827_REAL64, 3.147835047554940_REAL64, 7.267414688426937_REAL64, & ! 56
    5.157728581335382_REAL64, 7.906448979124780_REAL64, 2.044925813589564_REAL64, 6.781060775597956_REAL64, & ! 60 
    0.524855426380236_REAL64, 8.011723623624988_REAL64, 6.785686388193138_REAL64, 9.460089450840622_REAL64, & ! 64 
    0.915581428748301_REAL64, 9.084383705781780_REAL64, 5.099530370532040_REAL64, 6.149035651400551_REAL64, & ! 68 
    3.160712008859598_REAL64, 0.774874815154217_REAL64, 8.506140660037067_REAL64, 1.445269069485120_REAL64, & ! 72 
    3.704857926540119_REAL64, 6.223914330812510_REAL64, 9.975519024313703_REAL64, 5.173441332456926_REAL64, & ! 76 
    9.905112151503154_REAL64, 2.265344625558386_REAL64, 3.980051864202140_REAL64, 6.965686895116399_REAL64, & ! 80 
    0.646407594535227_REAL64, 7.476615823783636_REAL64, 4.204004312328131_REAL64, 8.113174202380232_REAL64, & ! 84
    3.796052857706114_REAL64, 3.190678202167843_REAL64, 9.860510667266144_REAL64, 7.181809052097272_REAL64, & ! 88 
    4.131833585698974_REAL64, 0.986302372429138_REAL64, 7.345591251664096_REAL64, 6.373062903534827_REAL64, & ! 92 
    0.738418767360661_REAL64, 1.205081664019784_REAL64, 9.815962089163827_REAL64, 4.967994217791492_REAL64]   ! 96

  ! Double Complex 
  zArray = [ &
    (2.599891163224722_REAL64, 7.868596456035510_REAL64), (5.115822433803345_REAL64, 5.625271533899060_REAL64), & ! 2
    (6.847937641435276_REAL64, 0.923965953476921_REAL64), (8.725789670394450_REAL64, 9.429378802197043_REAL64), & ! 4
    (0.965942366019155_REAL64, 8.458825052214273_REAL64), (9.093956148560236_REAL64, 0.113411405597612_REAL64), & ! 6
    (5.236827037454540_REAL64, 6.503452736589439_REAL64), (3.851453507059333_REAL64, 6.493020494573898_REAL64), & ! 8
    (7.628534797523447_REAL64, 5.756859813831817_REAL64), (6.319191487122334_REAL64, 2.782040038870991_REAL64), & ! 10
    (8.398364227998110_REAL64, 4.268349562390741_REAL64), (6.316223783496630_REAL64, 8.334669815207677_REAL64), & ! 12
    (2.701854636145602_REAL64, 4.008007361722771_REAL64), (5.542545924136538_REAL64, 4.438653568326196_REAL64), & ! 14
    (0.903841220286110_REAL64, 7.443813099870736_REAL64), (0.326154493029029_REAL64, 4.297425434330213_REAL64), & ! 16
    (0.372622794740213_REAL64, 9.757977530361101_REAL64), (5.223399080130893_REAL64, 9.096299785978145_REAL64), & ! 18
    (3.832479213400146_REAL64, 8.844522369047823_REAL64), (2.550176093311532_REAL64, 9.090498419428014_REAL64), & ! 20
    (8.945604447512554_REAL64, 3.985171609905782_REAL64), (6.250200018360077_REAL64, 5.675972861858413_REAL64), & ! 22
    (8.945118854449973_REAL64, 2.141658038678939_REAL64), (0.038591588607206_REAL64, 8.805814610317103_REAL64), & ! 24
    (2.351216350405121_REAL64, 2.448646223470780_REAL64), (6.409170351397391_REAL64, 3.045244812084679_REAL64), & ! 26
    (8.256208118888132_REAL64, 8.836887109472608_REAL64), (9.453729579854587_REAL64, 3.907953384046205_REAL64), & ! 28
    (8.013200647116648_REAL64, 1.571132879774678_REAL64), (6.251653967799038_REAL64, 6.989853527596546_REAL64), & ! 30
    (0.858689331204349_REAL64, 5.311795727884400_REAL64), (8.885665221405667_REAL64, 2.636660081896261_REAL64), & ! 32
    (2.347855253043302_REAL64, 8.396572277242370_REAL64), (4.955400297259374_REAL64, 1.523663664129055_REAL64), & ! 34
    (2.307681558094964_REAL64, 6.579539570846978_REAL64), (5.629483812400755_REAL64, 2.918287376488475_REAL64), & ! 36
    (6.223049847151239_REAL64, 7.159054082923781_REAL64), (2.807310340801989_REAL64, 4.122731585490750_REAL64), & ! 38
    (3.622057134141290_REAL64, 7.813920220025286_REAL64), (1.354871513019811_REAL64, 9.020695145496695_REAL64), & ! 40
    (2.896360858781397_REAL64, 4.995579261098015_REAL64), (7.835903881780232_REAL64, 6.770620626976057_REAL64), & ! 42
    (1.498117300188997_REAL64, 6.966194796179114_REAL64), (1.290130431227449_REAL64, 9.459462874157145_REAL64), & ! 44
    (8.864118492428657_REAL64, 5.150046129374196_REAL64), (6.794087459316357_REAL64, 9.767933578336322_REAL64), & ! 46
    (1.254603566400302_REAL64, 7.522430512464036_REAL64), (8.270519094561582_REAL64, 7.814303269102440_REAL64), & ! 48
    (1.908826667582241_REAL64, 4.286413406039776_REAL64), (0.144567471162916_REAL64, 3.252854887196962_REAL64), & ! 50
    (1.347020286263703_REAL64, 4.505174620840213_REAL64), (5.722749301785769_REAL64, 7.920232819831775_REAL64), & ! 52
    (4.197366099780537_REAL64, 5.325365443409556_REAL64), (9.257043920280161_REAL64, 8.990818377421817_REAL64), & ! 54
    (5.448333984687419_REAL64, 9.011239501703685_REAL64), (0.518267089040074_REAL64, 8.086025186801070_REAL64), & ! 56
    (3.349045793776864_REAL64, 2.286831172486750_REAL64), (8.224018076457371_REAL64, 3.482345789934266_REAL64), & ! 58
    (1.654705953682333_REAL64, 0.281335680082572_REAL64), (9.553697652775563_REAL64, 6.802907117056004_REAL64), & ! 60
    (8.605621147246641_REAL64, 9.390935153524339_REAL64), (6.801939644480517_REAL64, 9.174237484474340_REAL64), & ! 62
    (2.566916709011578_REAL64, 8.856181344290192_REAL64), (9.200427901181818_REAL64, 3.000634400366826_REAL64), & ! 64
    (0.733909466937926_REAL64, 7.673987349601374_REAL64), (0.849522158650140_REAL64, 7.287640795926718_REAL64), & ! 66
    (4.478896435877100_REAL64, 6.512486632169324_REAL64), (1.695018706149451_REAL64, 5.314492413595802_REAL64), & ! 68
    (6.338008581939598_REAL64, 0.140960041826924_REAL64), (4.703713352310025_REAL64, 8.863259711681058_REAL64), & ! 70
    (1.140284356780709_REAL64, 4.425411516501793_REAL64), (6.595478881615526_REAL64, 2.947729702676847_REAL64), & ! 72
    (9.503684015148998_REAL64, 6.942862475655555_REAL64), (2.068065270209161_REAL64, 5.547623259567388_REAL64), & ! 74
    (8.792784845119186_REAL64, 5.578576907197959_REAL64), (7.523340601514114_REAL64, 8.949007922295081_REAL64), & ! 76
    (8.418441730648903_REAL64, 1.308566826459481_REAL64), (1.891539338519886_REAL64, 1.536398308925564_REAL64), & ! 78
    (0.289020656650996_REAL64, 0.090848984918137_REAL64), (5.964517419720373_REAL64, 6.090493243838076_REAL64), & ! 80
    (9.189234320105523_REAL64, 7.335744433731223_REAL64), (3.011469521008484_REAL64, 4.955744790304295_REAL64), & ! 82
    (2.581630234983811_REAL64, 7.328543638765527_REAL64), (1.167607510415154_REAL64, 7.460415632576016_REAL64), & ! 84
    (8.097891070941914_REAL64, 7.452336812116621_REAL64), (3.371433141133766_REAL64, 5.843249804363376_REAL64), & ! 86
    (4.689521049559393_REAL64, 0.872650344600689_REAL64), (8.287174261600724_REAL64, 6.859447084388446_REAL64), & ! 88
    (2.673250550957546_REAL64, 9.694835645924528_REAL64), (1.837742026658353_REAL64, 2.999409229448979_REAL64), & ! 90
    (4.111848906888626_REAL64, 2.364891871154207_REAL64), (1.950553226810852_REAL64, 7.053792796459152_REAL64), & ! 92
    (1.805482536658427_REAL64, 5.223329179218547_REAL64), (2.961718232860060_REAL64, 4.627818882843439_REAL64), & ! 94
    (9.252324166792809_REAL64, 2.158891422960256_REAL64), (0.010088785070125_REAL64, 9.066062431168417_REAL64)]   ! 96

!=====================================================================
! Convert Arrays to Matrices
!=====================================================================
                    
   do colIndx=1, numColumns 
     do rowIndx=1, numRows
        linIndx = rowIndx + numRows*(colIndx-1)
        ! Single Real
        sMatrix(rowIndx, colIndx) = sArray(linIndx)
        ! Double Real
        dMatrix(rowIndx, colIndx) = dArray(linIndx)
        ! Single Complex
        cMatrix(rowIndx, colIndx) = cArray(linIndx) 
        ! Double Complex
        zMatrix(rowIndx, colIndx) = zArray(linIndx)
     enddo
   enddo 

!=====================================================================
! Make Outputs Doubly Singular (Mathematically Two Zeros)
!=====================================================================

  sArraySingular = sArray
  dArraySingular = dArray
  cArraySingular = cArray
  zArraySingular = zArray

  sMatrixSingular = sMatrix
  dMatrixSingular = dMatrix
  cMatrixSingular = cMatrix
  zMatrixSingular = zMatrix

  ! sArraySingular - Col 1 = Col 4 = Col 6
  sArraySingular(37:48) = sArraySingular(1:12)
  sArraySingular(61:72) = sArraySingular(1:12)

  ! dArraySingular - Col 5 = Col 7 = Col 8
  dArraySingular(73:84) = dArraySingular(49:60)
  dArraySingular(85:96) = dArraySingular(49:60)

  ! cArraySingular - Col 1 = Col 7 = Col 8
  cArraySingular(73:84) = cArraySingular(1:12)
  cArraySingular(85:96) = cArraySingular(1:12)

  ! zArraySingular - Col 3 = Col 4 = Col 7
  zArraySingular(37:48) = zArraySingular(25:36)
  zArraySingular(73:84) = zArraySingular(25:36)

  ! sMatrixSingular - Col 2 = Col 4 = Col 8
  sMatrixSingular(:, 4) = sMatrixSingular(:, 2)
  sMatrixSingular(:,  = sMatrixSingular(:, 2)

  ! dMatrixSingular - Col 2 = Col 5 = Col 7
  dMatrixSingular(:, 5) = dMatrixSingular(:, 2)
  dMatrixSingular(:, 7) = dMatrixSingular(:, 2)

  ! cMatrixSingular - Col 3 = Col 4 = Col 5
  cMatrixSingular(:, 4) = cMatrixSingular(:, 3)
  cMatrixSingular(:, 5) = cMatrixSingular(:, 3)

  ! zMatrixSingular - Col 1 = Col 3 = Col 8
  zMatrixSingular(:, 3) = zMatrixSingular(:, 1)
  zMatrixSingular(:,  = zMatrixSingular(:, 1)

!=====================================================================
! Query LAPACK Version
!=====================================================================

  call ilaver(LMajorVer, LMinorVer, LPatchVer)

!=====================================================================
! Check for Doubly Singular LUs
!=====================================================================

  write (*,*) '====================================================================='
  write (*, '(A,I0,A,I0,A,I0,A)') ' STARS_getrf_CheckSingular : Slightly Singular LUs (LAPACK ', &
            LMajorVer, '.', LMinorVer, '.', LPatchVer, ')'
  write (*,*) '====================================================================='

  warningStatus = 0

  ! Single Real Array
  call sgetrf(numRows, numColumns, sArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single real array doubly singular U     : Singular detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Single real array doubly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                              sArraySingular((linIndx-1)*12 + linIndx)
    enddo
  endif

  call sgetrf(numRows, numColumns, sMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single real matrix doubly singular U    : Singular Detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Single real matrix doubly singular U    : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16)') '     U(', linIndx, ', ', linIndx, ') = ', sMatrixSingular(linIndx, linIndx)
    enddo
  endif

  ! Double Real
  call dgetrf(numRows, numColumns, dArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double real array doubly singular U     : Singular detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Double real array doubly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                              dArraySingular((linIndx-1)*12 + linIndx)
    enddo
  endif

  call dgetrf(numRows, numColumns, dMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double real matrix doubly singular U    : Singular Detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Double real matrix doubly singular U    : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24)') '     U(', linIndx, ', ', linIndx, ') = ', dMatrixSingular(linIndx, linIndx)
    enddo
  endif

  ! Single Complex
  call cgetrf(numRows, numColumns, cArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single complex array doubly singular U  : Singular detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Single complex array doubly singular U  : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16,F21.16,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(cArraySingular((linIndx-1)*12 + linIndx)), &
                                                 AIMAG(cArraySingular((linIndx-1)*12 + linIndx)), 'i'
    enddo
  endif

  call cgetrf(numRows, numColumns, cMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single complex matrix doubly singular U : Singular Detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Single complex matrix doubly singular U : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16,F21.16,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(cMatrixSingular(linIndx, linIndx)), &
                                                 AIMAG(cMatrixSingular(linIndx, linIndx)), 'i'
    enddo
  endif

  ! Double Complex
  call zgetrf(numRows, numColumns, zArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double complex array doubly singular U  : Singular detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Double complex array doubly singular U  : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24,F29.24,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(zArraySingular((linIndx-1)*12 + linIndx)), &
                                                 AIMAG(zArraySingular((linIndx-1)*12 + linIndx)), 'i'
    enddo
  endif

  call zgetrf(numRows, numColumns, zMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double complex matrix doubly singular U : Singular Detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Double complex matrix doubly singular U : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24,F29.24,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(zMatrixSingular(linIndx, linIndx)), &
                                                 AIMAG(zMatrixSingular(linIndx, linIndx)), 'i'
    enddo
  endif

!=====================================================================
! Make Outputs Mostly Singular
!=====================================================================

  sArraySingular = sArray
  dArraySingular = dArray
  cArraySingular = cArray
  zArraySingular = zArray

  sMatrixSingular = sMatrix
  dMatrixSingular = dMatrix
  cMatrixSingular = cMatrix
  zMatrixSingular = zMatrix

  do colIndx=2, numColumns
    sArraySingular( (colIndx-1)*numRows + 1 : (colIndx-1)*numRows + 12 ) = sArray(1:12)
    dArraySingular( (colIndx-1)*numRows + 1 : (colIndx-1)*numRows + 12 ) = dArray(1:12)
    cArraySingular( (colIndx-1)*numRows + 1 : (colIndx-1)*numRows + 12 ) = cArray(1:12)
    zArraySingular( (colIndx-1)*numRows + 1 : (colIndx-1)*numRows + 12 ) = zArray(1:12)

    sMatrixSingular(:, colIndx) = sMatrix(:, 1)
    dMatrixSingular(:, colIndx) = dMatrix(:, 1)
    cMatrixSingular(:, colIndx) = cMatrix(:, 1)
    zMatrixSingular(:, colIndx) = zMatrix(:, 1)
  enddo

!=====================================================================
! Check for Mostly Singular LUs
!=====================================================================

  write (*,*) '====================================================================='
  write (*, ' (A,I0,A,I0,A,I0,A)') ' STARS_getrf_CheckSingular : Mostly Singular LUs (LAPACK ', &
            LMajorVer, '.', LMinorVer, '.', LPatchVer, ')'
  write (*,*) '====================================================================='
    
  errorStatus = 0

  ! Single Real
  call sgetrf(numRows, numColumns, sArraySingular, LDA, IPIV, INFO)


  if (INFO.ne.0) then
    write (*, *) 'Single real array mostly singular U     : Singular detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Single real array mostly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16)') '     U(', linIndx, ', ', linIndx, ') = ', &
                sArraySingular((linIndx-1)*12 + linIndx)
    enddo 
  endif

  call sgetrf(numRows, numColumns, sMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single real matrix mostly singular U    : Singular Detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Single real matrix mostly singular U    : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16)') '     U(', linIndx, ', ', linIndx, ') = ', sMatrixSingular(linIndx, linIndx)
    enddo 
  endif

  ! Double Real
  call dgetrf(numRows, numColumns, dArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double real array mostly singular U     : Singular detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Double real array mostly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                              dArraySingular((linIndx-1)*12 + linIndx)
    enddo
  endif

  call dgetrf(numRows, numColumns, dMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double real matrix mostly singular U    : Singular Detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Double real matrix mostly singular U    : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24)') '     U(', linIndx, ', ', linIndx, ') = ', dMatrixSingular(linIndx, linIndx)
    enddo
  endif

  ! Single Complex

  call cgetrf(numRows, numColumns, cArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single complex array mostly singular U  : Singular detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Single complex array mostly singular U  : Singular not detected ; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16,F21.16,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(cArraySingular((linIndx-1)*12 + linIndx)), &
                                                 AIMAG(cArraySingular((linIndx-1)*12 + linIndx)), 'i'
    enddo
  endif

  call cgetrf(numRows, numColumns, cMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single complex matrix mostly singular U : Singular Detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Single complex matrix mostly singular U : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16,F21.16,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(cMatrixSingular(linIndx, linIndx)), &
                                                 AIMAG(cMatrixSingular(linIndx, linIndx)), 'i'
    enddo
  endif

  ! Double Complex
  call zgetrf(numRows, numColumns, zArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double complex array mostly singular U  : Singular detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Double complex array mostly singular U  : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24,F29.24,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(zArraySingular((linIndx-1)*12 + linIndx)), &
                                                 AIMAG(zArraySingular((linIndx-1)*12 + linIndx)), 'i'
    enddo
  endif

  call zgetrf(numRows, numColumns, zMatrixSingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double complex matrix mostly singular U : Singular Detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Double complex matrix mostly singular U : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24,F29.24,A)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                                 REAL(zMatrixSingular(linIndx, linIndx)), &
                                                 AIMAG(zMatrixSingular(linIndx, linIndx)), 'i'
    enddo
  endif

!=====================================================================
! Report Summary and Return
!=====================================================================

  write (*,*) '====================================================================='

  if (warningStatus.eq.0) then
    write (*, '(A)')        ' mkl_check_singular         : No warnings found.'
  else
    write(*, '(A, I1, A)') ' mkl_check_singular WARNING: ', warningStatus, ' of 8 slightly singular LUs missed.'
  endif

  if (errorStatus.eq.0) then
    write (*, '(A)')        ' mkl_check_singular         : No errors found.'
  else
    write (*, '(A, I1, A)') ' mkl_check_singular ERROR   : ', errorStatus, ' of 8 mostly singular LUs missed.'
  endif

  write (*,*) '====================================================================='

end program mkl_check_singular

!=====================================================================
! EOF
!=====================================================================

 

 

 

 

 




getrfSingularMWE.zip
0 Kudos

Link Copied

10 Replies
mecej4
Honored Contributor III
‎05-15-2024 06:23 AM
3,087 Views

Lines 34 to 37 and 237 have syntax errors (check the subscript expressions!), so your example code cannot be compiled and run. To avoid errors, zip the source file(s) and attach the zip file to your post.

It would facilitate the diagnosis of the error if your example restricted itself to one case of a singular matrix that is passed to SGETRF and the returned value from Lapack/MKL of INFO=0. It is distracting to see two different precisions, real/complex, 1-D array and 2-D array representations of matrices, etc.

0 Kudos
Copy link
Kevin_T_
Beginner
‎05-15-2024 03:46 PM
3,027 Views

My apologies, it looks like the forum software is trying to convert things into emojis.

I've attached a zip file to the original post containing both the original example and one pruned smaller to just focus on sgetrf and degtrf (smaller example also added to this post).

As an additional experiment, I have Machine A where things work and Machine X where things don't. I compiled using statically-linked libraries on both machines as:

gfortran getrfSingularMWE_2.f90 -I"${MKLROOT}/include" -m64  -Wl,--start-group ${MKLROOT}/lib/intel64/libmkl_gf_lp64.a ${MKLROOT}/lib/intel64/libmkl_gnu_thread.a ${MKLROOT}/lib/intel64/libmkl_core.a -Wl,--end-group -lgomp -lpthread -lm -ldl -o getrfSingularMWE_2.x


Both executables detect singular matrices on Machine A and fail to detect singular matrices on Machine X.

program mkl_check_singular

  use iso_fortran_env, only: REAL32, REAL64

  implicit none

!=====================================================================
! This program checks the reporting of singular LU'd matrices using all
! combinations of single/double, real/complex, and array/matrix for an
! 8x12 matrix.
!
! Singular matrices are created by assigning multiple columns to be the
! same (matrix(:, 2) = matrix(:, 1), etc.).
!
! Different LU algorithms can fail to find slightly singular matrices 
! (i.e. only a few identical columns) in different ways, failing to 
! recognize some while recognizing others that alternative algorithms 
! don't find. If an algorithm fails to recognize a singular U when all 
! columns are identical, however, that algorithm should be treated with 
! suspicion. 
!=====================================================================

  integer :: numRows = 12, numColumns = 8, LDA = 12, INFO, &
             linIndx, rowIndx, colIndx, LMajorVer, LMinorVer, LPatchVer, &
             warningStatus, errorStatus 

  integer, dimension(8) :: IPIV
        
  real(kind=REAL32),       dimension(96) :: sArray, sArraySingular
  real(kind=REAL64),       dimension(96) :: dArray, dArraySingular

!=====================================================================
! Setup Arrays
!=====================================================================

  ! Single Real
  sArray = [ &
    8.984861_REAL32, 1.181552_REAL32, 9.884179_REAL32, 5.399821_REAL32, 7.069174_REAL32, 9.994916_REAL32, & ! 6
    2.878493_REAL32, 4.145225_REAL32, 4.648399_REAL32, 7.639571_REAL32, 8.182040_REAL32, 1.002215_REAL32, & ! 12
    1.781170_REAL32, 3.596349_REAL32, 0.567047_REAL32, 5.218857_REAL32, 3.358490_REAL32, 1.756690_REAL32, & ! 18
    2.089467_REAL32, 9.051536_REAL32, 6.753912_REAL32, 4.684682_REAL32, 9.121325_REAL32, 1.040116_REAL32, & ! 24
    7.455461_REAL32, 7.362675_REAL32, 5.618614_REAL32, 1.841941_REAL32, 5.972114_REAL32, 2.999370_REAL32, & ! 30
    1.341229_REAL32, 2.126015_REAL32, 8.949417_REAL32, 0.714528_REAL32, 2.424866_REAL32, 0.537544_REAL32, & ! 36
    4.417221_REAL32, 0.132832_REAL32, 8.971914_REAL32, 1.966582_REAL32, 0.933705_REAL32, 3.073669_REAL32, & ! 42
    4.560577_REAL32, 1.016694_REAL32, 9.953897_REAL32, 3.320928_REAL32, 2.973468_REAL32, 0.620452_REAL32, & ! 48
    2.982440_REAL32, 0.463513_REAL32, 5.054281_REAL32, 7.614259_REAL32, 6.310700_REAL32, 0.898917_REAL32, & ! 54
    0.808624_REAL32, 7.772405_REAL32, 9.051347_REAL32, 5.337720_REAL32, 1.091542_REAL32, 8.258089_REAL32, & ! 60
    3.380977_REAL32, 2.939731_REAL32, 7.463134_REAL32, 0.103366_REAL32, 0.484473_REAL32, 6.679161_REAL32, & ! 66
    6.034680_REAL32, 5.261025_REAL32, 7.297094_REAL32, 7.072535_REAL32, 7.813771_REAL32, 2.879770_REAL32, & ! 72
    6.925320_REAL32, 5.566698_REAL32, 3.965208_REAL32, 0.615907_REAL32, 7.801755_REAL32, 3.375839_REAL32, & ! 78
    6.078659_REAL32, 7.412540_REAL32, 1.048132_REAL32, 1.278884_REAL32, 5.495401_REAL32, 4.852294_REAL32, & ! 84
    8.904757_REAL32, 7.989603_REAL32, 7.343411_REAL32, 0.513319_REAL32, 0.728853_REAL32, 0.885275_REAL32, & ! 90
    7.983509_REAL32, 9.430081_REAL32, 6.837156_REAL32, 1.320830_REAL32, 7.227245_REAL32, 1.103535_REAL32]   ! 96

  ! Double Real
  dArray = [ &
    9.053306072443640_REAL64, 6.401942483322220_REAL64, 1.629395714406698_REAL64, 5.659105910016573_REAL64, & ! 4 
    9.316161369771947_REAL64, 7.831018836846311_REAL64, 6.856871369789888_REAL64, 4.662192849827422_REAL64, & ! 8
    2.603180956665296_REAL64, 5.692681710770191_REAL64, 2.487706628880768_REAL64, 3.193015915208970_REAL64, & ! 12 
    9.108022322169406_REAL64, 8.852200409944105_REAL64, 7.945894241658521_REAL64, 9.258098071668771_REAL64, & ! 16 
    1.788399507758780_REAL64, 5.175411482364652_REAL64, 6.270053164093095_REAL64, 9.131823459009045_REAL64, & ! 20 
    6.639682538388079_REAL64, 3.891928229926131_REAL64, 7.400075814161521_REAL64, 8.176348597087308_REAL64, & ! 24 
    6.003447931764931_REAL64, 0.849970645027028_REAL64, 9.223580010900541_REAL64, 0.535978153165582_REAL64, & ! 28
    5.270249526773659_REAL64, 1.188532766197476_REAL64, 3.801429902410049_REAL64, 8.128325330151092_REAL64, & ! 32 
    2.440959155716900_REAL64, 8.844226590935749_REAL64, 7.126467947591497_REAL64, 3.781484380821433_REAL64, & ! 36 
    2.489196109959332_REAL64, 2.528537452392730_REAL64, 7.672435945269350_REAL64, 0.498618850472403_REAL64, & ! 40 
    6.852885676089091_REAL64, 6.202781221510013_REAL64, 7.466846299687251_REAL64, 9.772556483456693_REAL64, & ! 44 
    3.839135126606154_REAL64, 2.602056023137224_REAL64, 8.774695360141239_REAL64, 8.060959835106466_REAL64, & ! 48 
    4.611210797796994_REAL64, 0.909616879242212_REAL64, 5.642688252731981_REAL64, 1.873828891486068_REAL64, & ! 52 
    5.316895858947326_REAL64, 3.550333217733827_REAL64, 3.147835047554940_REAL64, 7.267414688426937_REAL64, & ! 56
    5.157728581335382_REAL64, 7.906448979124780_REAL64, 2.044925813589564_REAL64, 6.781060775597956_REAL64, & ! 60 
    0.524855426380236_REAL64, 8.011723623624988_REAL64, 6.785686388193138_REAL64, 9.460089450840622_REAL64, & ! 64 
    0.915581428748301_REAL64, 9.084383705781780_REAL64, 5.099530370532040_REAL64, 6.149035651400551_REAL64, & ! 68 
    3.160712008859598_REAL64, 0.774874815154217_REAL64, 8.506140660037067_REAL64, 1.445269069485120_REAL64, & ! 72 
    3.704857926540119_REAL64, 6.223914330812510_REAL64, 9.975519024313703_REAL64, 5.173441332456926_REAL64, & ! 76 
    9.905112151503154_REAL64, 2.265344625558386_REAL64, 3.980051864202140_REAL64, 6.965686895116399_REAL64, & ! 80 
    0.646407594535227_REAL64, 7.476615823783636_REAL64, 4.204004312328131_REAL64, 8.113174202380232_REAL64, & ! 84
    3.796052857706114_REAL64, 3.190678202167843_REAL64, 9.860510667266144_REAL64, 7.181809052097272_REAL64, & ! 88 
    4.131833585698974_REAL64, 0.986302372429138_REAL64, 7.345591251664096_REAL64, 6.373062903534827_REAL64, & ! 92 
    0.738418767360661_REAL64, 1.205081664019784_REAL64, 9.815962089163827_REAL64, 4.967994217791492_REAL64]   ! 96


!=====================================================================
! Make Outputs Doubly Singular (Mathematically Two Zeros)
!=====================================================================

  sArraySingular = sArray
  dArraySingular = dArray

  ! sArraySingular - Col 1 = Col 4 = Col 6
  sArraySingular(37:48) = sArraySingular(1:12)
  sArraySingular(61:72) = sArraySingular(1:12)

  ! dArraySingular - Col 5 = Col 7 = Col 8
  dArraySingular(73:84) = dArraySingular(49:60)
  dArraySingular(85:96) = dArraySingular(49:60)

!=====================================================================
! Query LAPACK Version
!=====================================================================

  call ilaver(LMajorVer, LMinorVer, LPatchVer)

!=====================================================================
! Check for Doubly Singular LUs
!=====================================================================

  write (*,*) '====================================================================='
  write (*, '(A,I0,A,I0,A,I0,A)') ' STARS_getrf_CheckSingular : Slightly Singular LUs (LAPACK ', &
            LMajorVer, '.', LMinorVer, '.', LPatchVer, ')'
  write (*,*) '====================================================================='

  warningStatus = 0

  ! Single Real Array
  call sgetrf(numRows, numColumns, sArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single real array doubly singular U     : Singular detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Single real array doubly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                              sArraySingular((linIndx-1)*12 + linIndx)
    enddo
  endif

  ! Double Real
  call dgetrf(numRows, numColumns, dArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double real array doubly singular U     : Singular detected.'
  else
    warningStatus = warningStatus + 1
    write (*, *) 'Double real array doubly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                              dArraySingular((linIndx-1)*12 + linIndx)
    enddo
  endif

!=====================================================================
! Make Outputs Mostly Singular
!=====================================================================

  sArraySingular = sArray
  dArraySingular = dArray

  do colIndx=2, numColumns
    sArraySingular( (colIndx-1)*numRows + 1 : (colIndx-1)*numRows + 12 ) = sArray(1:12)
    dArraySingular( (colIndx-1)*numRows + 1 : (colIndx-1)*numRows + 12 ) = dArray(1:12)
  enddo

!=====================================================================
! Check for Mostly Singular LUs
!=====================================================================

  write (*,*) '====================================================================='
  write (*, ' (A,I0,A,I0,A,I0,A)') ' STARS_getrf_CheckSingular : Mostly Singular LUs (LAPACK ', &
            LMajorVer, '.', LMinorVer, '.', LPatchVer, ')'
  write (*,*) '====================================================================='
    
  errorStatus = 0

  ! Single Real
  call sgetrf(numRows, numColumns, sArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Single real array mostly singular U     : Singular detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Single real array mostly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F20.16)') '     U(', linIndx, ', ', linIndx, ') = ', &
                sArraySingular((linIndx-1)*12 + linIndx)
    enddo 
  endif

  ! Double Real
  call dgetrf(numRows, numColumns, dArraySingular, LDA, IPIV, INFO)

  if (INFO.ne.0) then
    write (*, *) 'Double real array mostly singular U     : Singular detected.'
  else
    errorStatus = errorStatus + 1
    write (*, *) 'Double real array mostly singular U     : Singular not detected; U(i, i) values are:'
    do linIndx=1, 8
      write (*, '(A,I1,A,I1,A,F28.24)') '     U(', linIndx, ', ', linIndx, ') = ', &
                                              dArraySingular((linIndx-1)*12 + linIndx)
    enddo
  endif

!=====================================================================
! Report Summary and Return
!=====================================================================

  write (*,*) '====================================================================='

  if (warningStatus.eq.0) then
    write (*, '(A)')        ' mkl_check_singular         : No warnings found.'
  else
    write(*, '(A, I1, A)') ' mkl_check_singular WARNING: ', warningStatus, ' of 2 slightly singular LUs missed.'
  endif

  if (errorStatus.eq.0) then
    write (*, '(A)')        ' mkl_check_singular         : No errors found.'
  else
    write (*, '(A, I1, A)') ' mkl_check_singular ERROR   : ', errorStatus, ' of 2 mostly singular LUs missed.'
  endif

  write (*,*) '====================================================================='

end program mkl_check_singular

!=====================================================================
! EOF
!=====================================================================






0 Kudos
Copy link
mecej4
Honored Contributor III
‎05-18-2024 05:31 AM
2,986 Views
In response to Kevin_T_

It would be helpful to have as detailed a description of the misbehaving Machine-X as possible -- CPU model, RAM size, OS, OneAPI and Compiler versions -- I imagine going to a PC store/merchant and asking if they can sell me a PC on which my test program will definitely fail!

On my Windows PC, with an AMD 4800U CPU, Win-11-64, OneAPI 2024.1, the shorter program ran and reported "Singular detected" in all four cases. Also with Gfortran, Cygwin, Lapack 3.11, on the same PC.

In response to Kevin_T_
0 Kudos
Copy link
Kevin_T_
Beginner
‎05-21-2024 05:04 PM
2,840 Views

Working Machine:

  1. i7-9750H, Linux Mint 19.1

Non-Working Machines:

  1. Xeon Gold 6226R, CentOS 7
  2. Xeon E5-2660 V4, CentOS 7
  3. Xeon Gold 6148, CentOS 7
  4. Xeon E5-2650 v3, Rocky Linux 8.9

I've tried the latest MKL standalone, 2024.1.0.695, as well as the MKL that ships with OneAPI 2022.2, 2023.2, and 2024.1. Also the MKL that shipped with Intel Parallel Studio 2018 Update 3. None of them seem to work properly on the affected machines.

I thought it might be a CentOS 7 thing, but if I read the documentation correctly OneAPI 2022 should still be supported, and I would have thought that at least Parallel Studio 2018 would be OK.

0 Kudos
Copy link
Kevin_T_
Beginner
‎05-30-2024 04:05 PM
530 Views

Any thoughts?

0 Kudos
Copy link
Ruqiu_C_Intel
Moderator
‎06-06-2024 11:42 PM
457 Views

Hi Kevin,


They run out " mkl_check_singular WARNING" and "mkl_check_singular ERROR" on my Xeon Icelake and core i7-12700 environment for both file getrfSingularMWE.f90 and getrfSingularMWE_2.f90 . My compiler and oneMKL version used in oneAPI 2024.1 tool kit package.


Here is my some logs:

getrfSingularMWE.f90 :

 =====================================================================

 mkl_check_singular WARNING: 7 of 8 slightly singular LUs missed.

 mkl_check_singular ERROR  : 2 of 8 mostly singular LUs missed.

 =====================================================================


getrfSingularMWE_2.f90:

 =====================================================================

 mkl_check_singular WARNING: 2 of 2 slightly singular LUs missed.

 mkl_check_singular ERROR  : 1 of 2 mostly singular LUs missed.

 =====================================================================


Can you share your whole logs from Working Machine and Non-Working Machines? And it's always good if you also share the logs when MKL_VERBOSE=1



Regards,

Ruqiu


0 Kudos
Copy link
Kevin_T_
Beginner
‎06-11-2024 03:12 PM
412 Views

Hi Ruqio,

I'm not sure what logs you might be referring to? At any rate an example output with "MKL_VERBOSE=1" is, on a working machine:


MKL_VERBOSE oneMKL 2023.0 Update 2 Product build 20230613 for Intel(R) 64 architecture Intel(R) Advanced Vector Extensions (Intel(R) AVX) enabled processors, Lnx 2.60GHz lp64 gnu_thread


and various combinations of ?getrf with:

MKL_VERBOSE DGETRF(12,8,0x7fff30198220,12,0x7fff301981f0,7) 1.79us CNR:OFF Dyn:1 FastMM:1 TID:0 NThr:4

An example on a non-working machine is:

MKL_VERBOSE oneMKL 2022.0 Update 1 Product build 20220311 for Intel(R) 64 architecture Intel(R) Advanced Vector Extensions 2 (Intel(R) AVX2) enabled processors, Lnx 2.00GHz lp64 gnu_thread

and various combinations of ?getrf with:

MKL_VERBOSE DGETRF(12,8,0x7ffff9ff8ff0,12,0x7ffff9ff89c0,0) 2.75us CNR:OFF Dyn:1 FastMM:1 TID:0 NThr:28

(Apologies for poor formatting, this forum's software is very broken in Edge on Windows 10).

 

0 Kudos
Copy link
Kevin_T_
Beginner
‎06-11-2024 03:23 PM
406 Views

I can confirm that exporting "MKL_ENABLE_INSTRUCTIONS=AVX" seems to fix the issue on the broken machines, implying that MKL's AVX2 code is unreliable.

0 Kudos
Copy link
Ruqiu_C_Intel
Moderator
‎06-13-2024 12:45 AM
343 Views

Hi Kevin,

 

Noted oneMKL version is different on your two machine. Working machine is using newer oneMKL version 2023.0 than non-working machine(oneMKL 2022.0).

Can you try latest oneMKL (version 2024.1) ?

Attached my log with AVX2 for new oneMKL version. 

 

Regards,

Ruqiu

 

 

 

 

 

 

 

getrfSingularMWE.log
getrfSingularMWE_2.log
0 Kudos
Copy link
Suheb
Beginner
‎06-13-2024 06:36 AM
331 Views

Oh no, that sounds really frustrating! Dealing with `getrf` routines not catching singular matrices can be a major headache, especially when it happens consistently across different setups. It’s good that you’ve got a minimal working example to illustrate the issue. Maybe there’s something subtle going on with the way the matrices are being handled or a quirk in the specific MKL versions you're using. Double-checking the example and ensuring all environment details are correct could help. Hopefully, someone with experience in this area can offer some insights or a workaround.

0 Kudos
Copy link
Reply

Community support is provided Monday to Friday. Other contact methods are available here.

Intel does not verify all solutions, including but not limited to any file transfers that may appear in this community. Accordingly, Intel disclaims all express and implied warranties, including without limitation, the implied warranties of merchantability, fitness for a particular purpose, and non-infringement, as well as any warranty arising from course of performance, course of dealing, or usage in trade.

For more complete information about compiler optimizations, see our Optimization Notice.