correctness (fm813 fm815 fm817 fm820 fm828 fm831) -- complex intrinsic

[schweitzpgi]

./fm831.out > fm831.log
Fortran STOP
IEEE arithmetic exceptions signaled: INEXACT
diff fm831.log .//ref/fm831.log
31,33c31
<     11     FAIL  
<                  COMPUTED= (     0.28402,     -0.99986)
<                  CORRECT=  (    -1.06820,     -0.03533)
---
>     11     PASS  
38,39c36,37
<                          11 TESTS PASSED
<                           1 TESTS FAILED
---
>                          12 TESTS PASSED
>                           0 TESTS FAILED

[sscalpone]

% cat c831.f
        COMPLEX AVC, CVC
        CVC = (3.125, 1.5)
        PRINT *, CVC
        PRINT *, COS(CVC)
        AVC = COS(CVC)
        PRINT *, AVC
      END

% flang c831.f
% a.out
 (3.125,1.5)
 (-.9998624,0.)
 (-.9998624,0.)

% gfortran c831.f
% a.out
             (3.12500000,1.50000000)
      (-2.35208583,-3.532877192E-02)
      (-2.35208583,-3.532877192E-02)

[jeanPerier]

That's a runtime ABI issue, the declaration of __fc_cos_1 is probably not correct. The result on X86 64 is that __fc_cos_1 is called with real part in $xmm0 and imaginary part in $xmm1, while __fc_cos_1 expects both parts in $xmm0. Hence, the runtime is only seeing the real part en returns cos((3.125, 0.)) = -.998624.

Declaring __fc_cos_1 as func @__fc_cos_1(f64) -> f64 and inserting cast (see below) in FIR solves the issue in X86 64, but I am not sure it is OK to do so for all targets.

Manual fix on X86 64:
Change fir.cos.z4.z4 to:

func @__fc_cos_1(f64) -> f64 attributes {fir.runtime}
func @fir.cos.z4.z4(%arg0: !fir.complex<4>) -> !fir.complex<4> attributes {fir.intrinsic, linkName = "internal"} {
  %0 = fir.alloca f64
  %1 = fir.convert %0 : (!fir.ref<f64>) -> !fir.ref<!fir.complex<4>>
  fir.store %arg0 to %1 : !fir.ref<!fir.complex<4>>
  %3 = fir.load %0 : !fir.ref<f64>
  %4 = call @__fc_cos_1(%3) : (f64) -> f64
  fir.store %4 to %0 : !fir.ref<f64>
  %5 = fir.load %1 : !fir.ref<!fir.complex<4>>
  return %5 : !fir.complex<4>
}

However, using __mth_i_ccos instead of __fc_cos_1 also solves the issue even if __mth_i_ccos is declared func @__mth_i_ccos(!fir.complex<4>) -> !fir.complex<4>, so maybe we should target this instead to avoid ABI hassles.

__mth_i_ccos:
https://github.com/flang-compiler/flang/blob/cb42a171b6462b2f3acd3e1841a56a1ae08291ef/runtime/libpgmath/lib/common/ccos.c#L10

__fc_cos_1 definition in libpgmath:
https://github.com/flang-compiler/flang/blob/cb42a171b6462b2f3acd3e1841a56a1ae08291ef/runtime/libpgmath/lib/common/mth_128defs_init.c#L4640

[jeanPerier]

Well, we will have to find a better way to deal with C99 Complex representation and lowering if we want to use it, LLVM sadly does not have an intrinsic type for it, hence the front-ends have to deal with the Complex ABI hiccups:

float _Complex __fc_cos_1(float _Complex); is translated by clang -S -emit-llvm into:

• On X86 64: declare dso_local <2 x float> @__fc_cos_1(<2 x float>)
• On X86 32 (-m32 option) i64 @__fc_cos_1({ float, float }* byval({ float, float }) align 4)
• On Arm 64 V8 declare { float, float } @__fc_cos_1([2 x float])

This was tested with the following C program. At least on X86 64, the program runs as expected.

#include <complex.h>
#include <stdio.h>

float _Complex __fc_cos_1(float _Complex);

int main(void) {
  float _Complex z1, z2;
  z1 = 3.125 + I*1.5;
  z2 = __fc_cos_1(z1);
  printf("(%f, %f)\n", creal(z2), cimag(z2));
}

That just show that handling on float _Complex in llvm environment is really not trivial, and that the middle-end emitting the LLVM IR has to adapt the generated IR depending on the target ABI.

Since there is currently no notion of target in FIR, it is hard to deal with this properly. For now, we can save the trouble by using the __mth_i_xxx runtime for complex, but we will have to face this with the following user code:

COMPLEX FUNCTION MY_CCOS(Z) BIND(C)
  COMPLEX, VALUE :: Z
END FUNCTION

That is guaranteed to be interoperable with the C function: float _Complex my_ccos(float _Complex) per Fortran 2018 Standard section 18.3.6 point 2 (3) and (4) + table 18.2.

[sscalpone]

MLIR has a complex type, doesn't it? Is it possible to extend MLIR to handle the intrinsics & lower them from MLIR to LLVM?

[jeanPerier]

MLIR has a complex type, doesn't it? Is it possible to extend MLIR to handle the intrinsics & lower them from MLIR to LLVM?

It does, unfortunately it looks like it also does not have the C99 _Complex linkage properties.

func @__fc_cos_1(%arg0 : complex<f32>) -> complex<f32> {
  return %arg0 : complex<f32>
}

through mlir-opt --convert-std-to-llvm and mlir-translate -mlir-to-llvmir (on X86 64):

define { float, float } @__fc_cos_1({ float, float } %0) !dbg !3 {
  ret { float, float } %0, !dbg !7
}

This is exactly what we also have with !fir.complex and does not work properly with C99 _Complex.

[schweitzpgi]

We need to be able to handle a target triple in our tools and add handlers in Tilikum to lower to the target-specific special sauce. Presently, we do simply let MLIR do its thing here. We should also look at their LLVM IR code gen. It may or may not have the right knobs for this.

[schweitzpgi]

#338 #339 #340 #341 #342 #343 could all be merged with this one. Just appear to be testing different complex intrinsic calls.